View Full Version : A bit OT: a discrete form of grad(divU).

okapota

04-27-2003, 11:33 AM

hello. i just read a paper about modelling clouds, and im implementing it's algorithm. in order to do that, i need the discrete form of W = grad(divU), where U is the velocity vector. so, in the paper, the x component of W is presented, and it is noted that the y and z are computed in a similiar way, so i sat and wrote the equations, and ended up with a slightly different solution, so maybe im wrong, or maybe there is a mistake in the paper. so if you can take a look at that-

my solution -

Wx = 1/4[Ux(x+2,y,z)+Ux(x-2,y,z)-2Ux(x,y,z)]

+ 1/4[Uy(x+1,y+1,z)+Uy(x-1,y-1,z)-Uy(x-1,y+1,z)-Uy(x+1,y-1,z)+Uz(x+1,y,z+1)+Uz(x-1,y,z-1)-Uz(x-1,y,z+1)-Uz(x+1,y,z-1)].

the papers solution, the second term is identical-

Wx = 1/2[Ux(x+1,y,z)+Ux(x-1,y,z)-2Ux(x,y,z)]

+ 1/4[Uy(x+1,y+1,z)+Uy(x-1,y-1,z)-Uy(x-1,y+1,z)-Uy(x+1,y-1,z)+Uz(x+1,y,z+1)+Uz(x-1,y,z-1)-Uz(x-1,y,z+1)-Uz(x+1,y,z-1)].

my question is who is right and what should be Wy,Wz?

i will appriciate any answer, thnx.

graphicsMan

04-27-2003, 02:15 PM

Is grad(divU) the divergence of U?

If so, then the grad dot U =

d(Ux)/dx + d(Uy)/dy + d(Uz)/dz

which ends up being estimated by

(1/2)*(Ux(x+1,y,z)-Ux(x-1,y,z)) etc... (by central differencing)

Brian

okapota

04-28-2003, 12:54 AM

its the gradient of the divergence of U.

i know, what you wrote is what used to start with.

graphicsMan

04-28-2003, 08:47 AM

Divergence is a scalar quantity... so the gradient should be a single-valued function. I don't get how you are getting Wx,Wy, and Wz...

Maybe I'm missing something?

bakery2k

04-28-2003, 09:55 AM

He wants the DIVERGENCE OF the gradient of U. This is a vector.

okapota

04-28-2003, 10:12 AM

bakery is right.

graphicsMan

04-28-2003, 08:23 PM

Okay, the answer I get goes:

1/4[Ux(x+1,y,z)+Ux(x-1,y,z)-2Ux(x,y,z)]...

This is because the 1/4 comes from 1/(deltaX^2) where your deltaX is 2. With your formula, you should get 1/16...

Brian

okapota

04-29-2003, 12:20 AM

ok, so why in the paper it is 1/2?

and what should y/z look like?

graphicsMan

05-01-2003, 09:51 PM

Originally posted by okapota:

ok, so why in the paper it is 1/2?

and what should y/z look like?

Good question! I don't know why the paper has 1/2. I believe that y/z should be analogous.

okapota

05-02-2003, 12:05 AM

by analogous you mean like that -

Wy = 1/2[Uy(x,y+1,z)+Uy(x,y-1,z)-2Uy(x,y,z)]

+ 1/4[Ux(x+1,y+1,z)+Ux(x-1,y-1,z)-Ux(x-1,y+1,z)-Ux(x+1,y-1,z)+Uz(x,y+1,z+1)+Uz(x,y-1,z-1)-Uz(x,y-1,z+1)-Uz(x,y+1,z-1)]/.

?

graphicsMan

05-02-2003, 10:59 AM

Yes... except that I think the 1/2 should be a 1/4?

Basically, you're just taking a different partial derivative, so it should be the same computation just with respect to a different variable.

The paper is incorrect. They copied their formula off Yanagita & Kaneko's paper, which had the incorrect formula.

PS You might find vorticity confinement and stronger incompressibility helpful for implementing clouds.

[This message has been edited by tie (edited 05-02-2003).]

okapota

05-03-2003, 12:09 AM

by saying stronger incompressibility you mean doing viscosity and pressure seperatly, each with a poission iterative solver?

okapota

05-04-2003, 10:48 AM

ok, nevermind. what about this, in a slide from gdc 2003 it is said that the decrete version of a poission equation in 2D is

P'i,j = 0.25*(Pi+1,j + Pi-1,j + Pi,j+1 + Pi,j-1 - Delta^2(DivU).

i understand this, but i wont to implement such a solver for a 3d grid, and i dont know how to derive the formula, anybody knows?

Powered by vBulletin® Version 4.2.3 Copyright © 2017 vBulletin Solutions, Inc. All rights reserved.