<html>
<head>
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
</head>
<body bgcolor="#FFFFFF" text="#000000">
Hi Neks !<br>
<br>
We are interested in studying the weakly non linear dynamics of a
perturbation. We are therefore using Nek in perturbation mode, and
we have implemented the usually omitted 2nd order non linear
perturbation convection terms in advabp() in perturb.f as follows :<br>
<br>
<pre> subroutine advabp
C
C Eulerian scheme, add convection term to forcing function
C at current time step.
C
include 'SIZE'
include 'INPUT'
include 'SOLN'
include 'MASS'
include 'TSTEP'
C
COMMON /SCRNS/ TA1 (LX1*LY1*LZ1*LELV)
$ , TA2 (LX1*LY1*LZ1*LELV)
$ , TA3 (LX1*LY1*LZ1*LELV)
$ , TB1 (LX1*LY1*LZ1*LELV)
$ , TB2 (LX1*LY1*LZ1*LELV)
$ , TB3 (LX1*LY1*LZ1*LELV)
C
ntot1 = nx1*ny1*nz1*nelv
ntot2 = nx2*ny2*nz2*nelv
c
if (if3d) then
call opcopy (tb1,tb2,tb3,vx,vy,vz) ! Save velocity
call opcopy (vx,vy,vz,vxp(1,jp),vyp(1,jp),vzp(1,jp)) ! U <-- dU
call convop (ta1,tb1) ! du.grad U
call convop (ta2,tb2)
call convop (ta3,tb3)
call opcopy (vx,vy,vz,tb1,tb2,tb3) ! Restore velocity
c
do i=1,ntot1
tmp = bm1(i,1,1,1)*vtrans(i,1,1,1,ifield)
bfxp(i,jp) = bfxp(i,jp)-tmp*ta1(i)
bfyp(i,jp) = bfyp(i,jp)-tmp*ta2(i)
bfzp(i,jp) = bfzp(i,jp)-tmp*ta3(i)
enddo
c
call convop (ta1,vxp(1,jp)) ! U.grad dU
call convop (ta2,vyp(1,jp))
call convop (ta3,vzp(1,jp))
c
do i=1,ntot1
tmp = bm1(i,1,1,1)*vtrans(i,1,1,1,ifield)
bfxp(i,jp) = bfxp(i,jp)-tmp*ta1(i)
bfyp(i,jp) = bfyp(i,jp)-tmp*ta2(i)
bfzp(i,jp) = bfzp(i,jp)-tmp*ta3(i)
enddo
<font color="#cc0000">
! ADD NON LINEAR TERM
call opcopy (tb1,tb2,tb3,vx,vy,vz)
call opcopy (vx,vy,vz,vxp(1,jp),vyp(1,jp),vzp(1,jp))
call convop (ta1,vxp(1,jp))
call convop (ta2,vyp(1,jp))
call convop (ta3,vzp(1,jp))
call opcopy (vx,vy,vz,tb1,tb2,tb3)
do i=1,ntot1
tmp = bm1(i,1,1,1)*vtrans(i,1,1,1,ifield)
bfxp(i,jp) = bfxp(i,jp)-tmp*ta1(i)
bfyp(i,jp) = bfyp(i,jp)-tmp*ta2(i)
bfzp(i,jp) = bfzp(i,jp)-tmp*ta3(i)
enddo</font>
c
else ! 2D
c
call opcopy (tb1,tb2,tb3,vx,vy,vz) ! Save velocity
call opcopy (vx,vy,vz,vxp(1,jp),vyp(1,jp),vzp(1,jp)) ! U <-- dU
call convop (ta1,tb1) ! du.grad U
call convop (ta2,tb2)
call opcopy (vx,vy,vz,tb1,tb2,tb3) ! Restore velocity
c
do i=1,ntot1
tmp = bm1(i,1,1,1)*vtrans(i,1,1,1,ifield)
bfxp(i,jp) = bfxp(i,jp)-tmp*ta1(i)
bfyp(i,jp) = bfyp(i,jp)-tmp*ta2(i)
enddo
c
call convop (ta1,vxp(1,jp)) ! U.grad dU
call convop (ta2,vyp(1,jp))
c
do i=1,ntot1
tmp = bm1(i,1,1,1)*vtrans(i,1,1,1,ifield)
bfxp(i,jp) = bfxp(i,jp)-tmp*ta1(i)
bfyp(i,jp) = bfyp(i,jp)-tmp*ta2(i)
enddo
<font color="#cc0000"> ! ADD NON LINEAR TERM
call opcopy (tb1,tb2,tb3,vx,vy,vz)
call opcopy (vx,vy,vz,vxp(1,jp),vyp(1,jp),vzp(1,jp))
call convop (ta1,vxp(1,jp))
call convop (ta2,vyp(1,jp))
call opcopy (vx,vy,vz,tb1,tb2,tb3)
do i=1,ntot1
tmp = bm1(i,1,1,1)*vtrans(i,1,1,1,ifield)
bfxp(i,jp) = bfxp(i,jp)-tmp*ta1(i)
bfyp(i,jp) = bfyp(i,jp)-tmp*ta2(i)
enddo</font>
c
endif
c
return
end</pre>
<br>
Is this the correct implementation ? Are there any other
modifications to be made ? An alternative would obviously be to add
the non linear perturbation convection term in USERF, would that be
a safer solution ?<br>
<br>
Any help would be greatly appreciated ! <br>
<br>
Best, <br>
<br>
Holly<br>
</body>
</html>