Forum Message

 

 

We have moved the forum to https://forum.ngsolve.org . This is an archived version of the topics until 05/05/23. All the topics were moved to the new forum and conversations can be continued there. This forum is just kept as legacy to not invalidate old links. If you want to continue a conversation just look for the topic in the new forum.

Notice

The forum is in read only mode.

JacobianInverse as specialcf

More
1 year 5 months ago #4535 by creativeworker
Is there a way to access the JacobianInverse as specialcf in python?
More
1 year 5 months ago - 1 year 5 months ago #4536 by mneunteufel
Hi,

you can invert the specialcf.JacobianMatrix to obtain the inverse
F_inv = Inv(specialcf.JacobianMatrix(3))

Best,
Michael
Last edit: 1 year 5 months ago by mneunteufel.
The following user(s) said Thank You: creativeworker
More
1 year 5 months ago - 1 year 5 months ago #4541 by creativeworker
Perfect. Is there a way to Diff a JacobianMatrix?
The goal is to implement the virtual work principle:  www.sciencedirect.com/science/article/pii/0304885382901883
Last edit: 1 year 5 months ago by creativeworker.
More
1 year 5 months ago #4543 by mneunteufel
It should be possible to interpolate the Jacobian (or some components of it) into an L2/VectorL2/MatrixValued(L2) GridFunction and then compute the gradient.
The following user(s) said Thank You: creativeworker
More
1 year 5 months ago #4544 by creativeworker
 
Thanks, the MatrixValued(L2(mesh)) helped me out. I attach the minimal example for other users!
Code:
[code]# -*- coding: utf-8 -*- from ngsolve import * from netgen.geom2d import SplineGeometry # geometry geo = SplineGeometry() geo.AddRectangle(p1 = (-10,-5), p2 = (10, 5), leftdomain = 1, rightdomain = 0, bc = "outer") geo.AddRectangle(p1 = (-5,-3), p2 = (5,3), leftdomain = 2, rightdomain = 1, bc = 'gamma') geo.SetMaterial(1, 'mat1') geo.SetMaterial(2, 'mat2') ngmesh = geo.GenerateMesh(maxh = 5) mesh = Mesh(ngmesh) #%% Jacobian G = specialcf.JacobianMatrix(2) print('Jacobian Dims: ', G.dims.NumPy()) G_inv = Inv(G)  print('Jacobian Inv Dims: ', G_inv.dims.NumPy()) G_det = Det(G)  print('Jacobian det dims: ', G_det.dim) G_det_inv = 1/G_det mesh.UnsetDeformation() gfu_G = GridFunction(MatrixValued(L2(mesh))) gfu_G.Set(G) gfu_G_det = GridFunction(L2(mesh)) gfu_G_det.Set(G_det) #% Finite displacement vdisplacement = 0.5 fes_structure = H1(mesh, order=1, dim=mesh.dim) gfu_struc = GridFunction(fes_structure) gfu_struc.Set((vdisplacement,0), definedon=mesh.Materials('mat2')) mesh.SetDeformation(gfu_struc) gfu_G2 = GridFunction(MatrixValued(L2(mesh))) gfu_G2.Set(G) gfu_G2_det = GridFunction(L2(mesh)) gfu_G2_det.Set(G_det) mesh.UnsetDeformation() Gdx = (gfu_G2-gfu_G)/vdisplacement G_det_dx = (gfu_G2_det-gfu_G_det)/vdisplacement gfu_Gdet_dx = GridFunction(L2(mesh)) gfu_Gdet_dx.Set(G_det_dx) #  vtk = VTKOutput(ma = mesh,coefs=[gfu_G_det, gfu_G2_det, gfu_Gdet_dx, gfu_G, Gdx],names=['G_det','G_det2','G_detdx','G','Gdx'],filename = "simple_rectange",subdivision=0, legacy=True) vtk.Do()
[/code]
Time to create page: 0.152 seconds