modified: src/pressure.f90

src/pressure.f90

@@ -32,7 +32,7 @@
 
         call mass2(ra,i,j,l)
         de=  (defect_c(i,j,l)+ra)/h2sum*urfp
-        dpa= de*(1._pr/dt+2._pr*nuf/((xf(i)-xf(i-1))**2+(yf(j)-yf(j-1))**2+(zf(l)-zf(l-1))**2))
+        dpa= de*rhof*(1._pr/dt+2._pr*nuf/((xf(i)-xf(i-1))**2+(yf(j)-yf(j-1))**2+(zf(l)-zf(l-1))**2))
         dex= de/(xf(i)-xf(i-1))
         dey= de/(yf(j)-yf(j-1))
         dez= de/(zf(l)-zf(l-1))