← Back to branch summary

~cbc-core/cbc.solve/main

« back to all changes in this revision

Viewing changes to cbc/swing/fsinewton/solver/jacobianforms.py

  • Committer: Gabriel Balaban
  • Date: 2012-08-27 08:09:11 UTC
  • Revision ID: gabrib@math.uio.no-20120827080911-si5n21ghq44w3asu
some bugs in Newton's method analytic problem. Investigation required

Show diffs side-by-side

added added

removed removed

Lines of Context:
cbc/swing/fsinewton/solver/jacobianforms.py
197
197
        
198
198
    for dxF in dxFlist:
199
199
        #DT  
200
 
        J_FM =  inner(v_F, rho*J(D_F)*tr(dot(grad(dD_F), inv(F(D_F))))*dotU)*dxF
201
 
        J_FM +=  inner(v_F, rho*J(D_F)*tr(dot(grad(dD_F), inv(F(D_F))))*dot(grad(U), dot(inv(F(D_F)), U - dotD_F)))*dxF
 
200
        J_FM =  inner(v_F, rho*J(D_F)*tr(dot(grad(dD_F), inv(F(D_F))))*dotU)
 
201
        J_FM +=  inner(v_F, rho*J(D_F)*tr(dot(grad(dD_F), inv(F(D_F))))*dot(grad(U), dot(inv(F(D_F)), U - dotD_F)))
202
202
        J_FM += -inner(v_F,rho*J(D_F)*dot((dot(grad(U), dot(inv(F(D_F)), \
203
 
                 dot(grad(dD_F), inv(F(D_F)))))), U - dotD_F ))*dxF
204
 
        J_FM += -inner(v_F, rho*J(D_F)*dot(grad(U), dot(inv(F(D_F)),dotdD_F)))*dxF
 
203
                 dot(grad(dD_F), inv(F(D_F)))))), U - dotD_F ))
 
204
        J_FM += -inner(v_F, rho*J(D_F)*dot(grad(U), dot(inv(F(D_F)),dotdD_F)))
205
205
 
206
206
        #SigmaF
207
 
        J_FM += inner(grad(v_F),dD_FSigmaF(D_F,dD_F,U,P,mu))*dxF
 
207
        J_FM += inner(grad(v_F),dD_FSigmaF(D_F,dD_F,U,P,mu))
208
208
 
209
209
        #Div U_F (incompressibility)
210
 
        J_FM +=  inner(q, div(J(D_F)*tr(dot(grad(dD_F), inv(F(D_F))))*dot(inv(F(D_F)), U)))*dxF
211
 
        J_FM += -inner(q, div(J(D_F)*dot(dot(inv(F(D_F)), grad(dD_F)), dot(inv(F(D_F)), U))))*dxF
 
210
        J_FM +=  inner(q, div(J(D_F)*tr(dot(grad(dD_F), inv(F(D_F))))*dot(inv(F(D_F)), U)))
 
211
        J_FM += -inner(q, div(J(D_F)*dot(dot(inv(F(D_F)), grad(dD_F)), dot(inv(F(D_F)), U))))
 
212
        if dFSIlist is None:
 
213
            J_FM = J_FM*dxF
 
214
        else:
 
215
            J_FM = J_FM('+')*dxF
212
216
        forms.append(J_FM)
213
217
    return sum(forms)
214
218