Hi Florian,
Answers below.
On 9/18/07, Florian De Vuyst <[EMAIL PROTECTED]> wrote:
>
> Hi Daniel,
>
> Thanks a lot for your answer, as you ask me, I give the elesticity
> equations.
> Imagine in 2D, U=(u,v) a two-dimensional spatial field, mu and lambda
> are the Lamé parameters
> (constant or possibly spatial functions). F=(f,g) is an external force
> field.
>
> Then the system of two equations is say :
>
> div.[ mu grad(U) + mu [grad(U)]^T + lambda I tr(grad(U)) ] + F = (0,0).
>
> As you see, both u and v appear in the two equations,
> Inside the div.[...] you will matrices.
>
> In your opinion, can it be possible to embed those equations in the Fipy
> framework ?
Yes. Here's how. First create a mesh then create the two field variables
mesh = ...
u = CellVariable(mesh=mesh)
v = CellVariable(mesh=mesh)
The equations will be,
Xhat = (1, 0)
Yhat = (0, 1)
u_eqn = DiffusionTerm(mu) \
+ mu * (u.getFaceGrad().dot(Xhat) * Xhat \
+ v.getFaceGrad().dot(Xhat) * Yhat).getDivergence() \
+ ConvectionTerm(lambda * Xhat) \
+ lambda * v.getGrad().dot(Yhat) \
+ f
v_eqn = DiffusionTerm(mu) \
+ mu * (u.getFaceGrad().dot(Yhat) * Xhat \
+ v.getFaceGrad().dot(Yhat) * Yhat).getDivergence() \
+ ConvectionTerm(lambda * Yhat) \
+ lambda * u.getGrad().dot(Xhat) \
+ g
Then iterate the equations
for step in range(steps):
u_eqn.solve(u)
v_eqn.solve(v)
Note the following about the above equations:
* I haven't tested them. Get back to me if there are problems.
* I'm assuming constant mu and lambda in space.
* The lambda term has been made implicit, could do a similar thing with the
transposed term.
* Fipy would have a neater script if we could define U = (f, g)
and have a vector
equation. This is not possible yet. I could do with that in my
work as well. We have
vector variables, just not vector equations. This requires
fixing boundary conditions
as well.
* Boundary conditions could get a little tricky with the way the
above equations
are defined. Check back with me if you have problems.
* You may need some relaxation to solve to equilibrium.
Hope this helps.
Cheers
--
Daniel Wheeler
