Dear Yiannis, Thanks for your answer. Unfortunately it gives me the same matrix which is not positive definite. Developing your notation leads to my much longer expression (the terms j=k canceled out)
Best Regards Louis On 27 Nov 2013, at 15:28, Ioannis Koufogiannis <[email protected]> wrote: > Dear Louis, > > some time ago I developed such a code and if I remember well it was working > properly. > The piece of code that I used for assembling the stiffness matrix, using the > C++ getfem functions, was: > > getfem::generic_assembly matrix_A; > ... > ... > std::stringstream sstr; > sstr << "t=comp(vGrad(#1).vGrad(#1));""M(#1,#1)+=sym((t(:,j,k,:,j,k) - > t(:,j,k,:,k,j) - t(:,k,j,:,j,k) + t(:,k,j,:,k,j))*0.5)"; > matrix_A.set(sstr.str()); > > I arrived to this result after some expansions and rearrangements of the > vectorial products on the paper. > At the moment, I am not sure how this differs from your development, but you > can try it. > It seems that in my calculations, compared to yours, I have all the > combinations of (j,k), 1<j,k<3 and then I half the final result. > > Best, > Yiannis > > On 27.11.2013 14:27, Louis Kovalevsky wrote: >> Dear Getfem++ users, >> >> I would like to use Nedelec finite element to calculate the mode of 3d >> electromagnetic cavity. >> I'm using the matlab interface. >> The importation of the mesh, definition of element and integration method is >> fine: >> m = gf_mesh('import', 'gmsh', 'cylinder.msh'); >> mf = gf_mesh_fem(m,3); >> gf_mesh_fem_set(mf,'fem',gf_fem('FEM_NEDELEC(3)')); >> mim=gfMeshIm(m,gf_integ('IM_TETRAHEDRON(8)')); >> >> To assemble the term \int E.E'dx I'm using: >> Ma=gf_asm('volumic','M(#1,#1)+=comp(vBase(#1).vBase(#1))(:,k,:,k)', mim,mf); >> >> I got the same result with "gf_asm('mass matrix', mim, mf);" ( as excepted) >> >> For the "stiffness" matrix \int curlE curlE' dx i'm using: >> K=gf_asm('volumic','M(#1,#1)+=comp(vGrad(#1).vGrad(#1))(:,2,3,:,2,3)+comp(vGrad(#1).vGrad(#1))(:,3,2,:,3,2)-comp(vGrad(#1).vGrad(#1))(:,2,3,:,3,2)-comp(vGrad(#1).vGrad(#1))(:,3,2,:,2,3)+comp(vGrad(#1).vGrad(#1))(:,3,1,:,3,1)+comp(vGrad(#1).vGrad(#1))(:,1,3,:,1,3)-comp(vGrad(#1).vGrad(#1))(:,3,1,:,1,3)-comp(vGrad(#1).vGrad(#1))(:,1,3,:,3,1)+comp(vGrad(#1).vGrad(#1))(:,1,2,:,1,2)+comp(vGrad(#1).vGrad(#1))(:,2,1,:,2,1)-comp(vGrad(#1).vGrad(#1))(:,1,2,:,2,1)-comp(vGrad(#1).vGrad(#1))(:,2,1,:,1,2)', >> mim,mf); >> >> My problem is that the obtained matrix K is not positive definite (and it >> should be...). >> >> Does anyone have a suggestion why the K matrix is not positive definite? >> >> Thank you in advance for the answer. >> >> Best Regards >> >> Louis Kovalevsky >> >> >> >> >> >> _______________________________________________ >> Getfem-users mailing list >> [email protected] >> https://mail.gna.org/listinfo/getfem-users > > _______________________________________________ > Getfem-users mailing list > [email protected] > https://mail.gna.org/listinfo/getfem-users
_______________________________________________ Getfem-users mailing list [email protected] https://mail.gna.org/listinfo/getfem-users
