[Nek5000-users] Impose/change velocity in userchk

[nek5000-users at lists.mcs.anl.gov]

Dear Neks,

I am using the MHD example for kinematic Dynamo.
In the example gpf, the force is added as
f=-Nabla^2 u - partial u/partial t, where u is
         ux   = -sys
         uy   = -cxc
         uz   = sxc + cys

I want to do it the straight forward approach i.e. directly imposing the
velocity in userchk and switching off the iteration Navier-Stokes equation.
For this I set the IFNAV to False
  F T F F F F F F F F F IFNAV & IFADVC (convection in P.S. fields)
  F F F T T T T T T T T T IFTMSH (IF mesh for this field is T mesh)

IC being the same as imposed velocity

and in userchk added the following

      if (ifield.eq.1) then ! velocity
         omega = 1.
      do i=1,nx1*ny1*nz1*nelv
x = xm1(i,1,1,1); y = ym1(i,1,1,1); z = zm1(i,1,1,1)

         st  = sin(omega*time)
         ct  = cos(omega*time)
c
         argxc = x + ct
         argyc = y + ct
         argxs = x + st
         argys = y + st
c
         visc = param(2)
         s    = 3./2.
         s    = sqrt(s)
c
         cxc = s*cos(argxc)
         cys = s*cos(argys)
         sxc = s*sin(argxc)
         sys = s*sin(argys)
c
         vx(i,1,1,1) = -sys
         vy(i,1,1,1) = -cxc
         vz(i,1,1,1) = sxc + cys
      enddo
      endif

      ibstep = iostep
      if (istep.gt.0.and.mod(istep,ibstep).eq.0)
     $   call outpost(bx,by,bz,pm,t,'   ') ! odd blah.fld/.f for B-field
c    $   call outpost(bx,by,bz,pm,t,'mag') !  magblah.fld/.f for B-field
      if (ifxyo.and.istep.gt.istep) ifxyo = .false.

 kinerg2 = 0.5*(glsc3(vx,bm1,vx,m)
     $        +glsc3(vy,bm1,vy,m)
     $        +glsc3(vz,bm1,vz,m))/volvm1


      n = nx1*ny1*nz1*nelfld(ifldmhd) ! Magnetic energy !?nelb
      energy0 = energy
      energy  = glsc3(bx,bm1,bx,n)
     $        + glsc3(by,bm1,by,n)
     $        + glsc3(bz,bm1,bz,n)
      energy  = 0.5*energy/volfld(ifldmhd)

      if (istep.gt.1) then
         ratio  = energy/energy0
         growth = 0.5*log(ratio)/dt ! 1/2 of magnetic energy growth rate
      if (nid.eq.0) write(6,1) istep,time,growth,energy,
     $ kinerg1,kinerg2
    1    format(i8,1p5e15.6,' growth')

         nstepag = param(120) ! p121 < t_avg < p120
         nstepsa = param(121)
         if (istep.lt.nstepag) call
     $      runtimeavg(Emga,growth,1,nstepsa,1,'gr_Em') ! avg. growth
      endif

When I calculate the energy of the system, it is decaying.
Is there anyway I can tamper with / change the velocity field array in
nek5000 from userchk?

Thanking you in advance.

Cheers,
Sandeep

PS: Attached are the .usr and .rea files
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20161121/741a3caa/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: gpf.usr
Type: application/octet-stream
Size: 12212 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20161121/741a3caa/attachment-0002.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: gpf.rea
Type: application/octet-stream
Size: 5891 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20161121/741a3caa/attachment-0003.obj>