Hi Neks, [changing the subject from "Integration Method"]
I have tried to implement calculation of the pressure gradient and adding it to the particle motion (ppiclfT subroutine). BUt I have found the dp/dx=0, dpdy=0. Clearly I am doing something wrong. (I assume that the pressure field is "pr" ?)
Can anyone help out here.
Here are the parts of my code that pertains to this:
[1] in *.usr:
! Pressure-gradient field
call rzero(dpdx,n)
call rzero(dpdy,n)
call gradm1(dpdx,dpdy,dpdz,pr)
call ppiclf_solve_InterpFieldUser(PPICLF_R_JDPDX,dpdx(1,1,1,1))
call ppiclf_solve_InterpFieldUser(PPICLF_R_JDPDY,dpdy(1,1,1,1))
if(ppiclf_npart .ne. 0) then
do i = 1,ppiclf_npart
dpx = ppiclf_rprop(PPICLF_R_JDPDX,i)
dpy = ppiclf_rprop(PPICLF_R_JDPDY,i)
[2] in subroutine ppiclf_user_SetYdot:
! Pressure-gradient
fdpdx = -ppiclf_rprop(PPICLF_R_JVOLP,i)*ppiclf_rprop(PPICLF_R_JDPDX,i)
fdpdy = -ppiclf_rprop(PPICLF_R_JVOLP,i)*ppiclf_rprop(PPICLF_R_JDPDY,i)
! set ydot for all PPICLF_SLN number of equations
ppiclf_ydot(PPICLF_JVX,i) = (fbx+fqsx+fdpdx+fcx)/rmass
ppiclf_ydot(PPICLF_JVY,i) = (fby+fqsy+fdpdy+fcy)/rmass
ppiclf_ydotc(PPICLF_JVX,i) = -fqsx-fdpdx
ppiclf_ydotc(PPICLF_JVY,i) = -fqsy-fdpdy
Thanks