From ac43d3e3878e7e207c1730011a35a2fd89d25ec5 Mon Sep 17 00:00:00 2001
From: Dave Whipp <dwhipp@dal.ca>
Date: Wed, 8 Jun 2011 10:28:09 +0000
Subject: [PATCH] First attempt at implementing the surface stabilization of
Kaus et al., 2010
---
src/build_system_wsmp.f90 | 15 ++++++++++++++-
1 file changed, 14 insertions(+), 1 deletion(-)
diff --git a/src/build_system_wsmp.f90 b/src/build_system_wsmp.f90
index 8ec16c1f..95dd444e 100644
--- a/src/build_system_wsmp.f90
+++ b/src/build_system_wsmp.f90
@@ -68,7 +68,7 @@ double precision weightel_glob(osolve%nleaves)
integer, external :: ifind_loc_wsmp
double precision,dimension(:,:),allocatable::lsf_el
-double precision ael(8*ndof,8*ndof),bel(8*ndof),patch1(5,5),patch2(3,3)
+double precision ael(8*ndof,8*ndof),bel(8*ndof),patch1(5,5),patch2(3,3),elem_dz
double precision,dimension(:),allocatable::penalty_loc,penalty_glob,eviscosity,weightel
double precision r0,s0,t0,rst,amin
integer i,j,k,itodo,idof,i1,i2,iii,jjj,err,jj,i1loc,k1,k2,ii
@@ -172,6 +172,8 @@ do ileaves=1,osolve%nleaves
! transfer elemental matrix and rhs into wsmpmatrix structure
! different for symmetric or nonsymmetric matrices
+ do i=1,params%mpe
+
do k1=1,params%mpe
do idof1=1,ndof
i1=(vo%ftr(osolve%icon(k1,ileaves))-1)*ndof+idof1
@@ -189,6 +191,17 @@ do ileaves=1,osolve%nleaves
endif
enddo
enddo
+ if (kfix(ileaves) < 1) then
+ if (idof1 == 3) then
+ iloc=iffind_loc_wsmp(i1loc,i1,ja,ia,n_iproc,nz_loc)
+ elem_dz=(osolve%icon(5,ileaves)-osolve%icon(1,ileaves))*params%vex
+ if (k1 < 5) then
+ avals(iloc)=avals(iloc)+(params%dt*weightel(ileaves))/(2*elem_dz)
+ else
+ avals(iloc)=avals(iloc)-(params%dt*weightel(ileaves))/(2*elem_dz)
+ endif
+ endif
+ endif
b(i1loc,1)=b(i1loc,1)+bel(iii)
endif
enddo
--
GitLab