|
From: | Louis Kovalevsky |
Subject: | Re: [Getfem-users] Assembling wave equation for calculating eigenmodes |
Date: | Sun, 1 Dec 2013 15:43:08 +0000 |
Dear Johnathan, When using finite element, you are solving the governing equation (in your case, the Maxwell equation) using a weak formulation. The discretisation of the unknown field injected into the weak formulation leads to an algebraical system. It seems to me that you're not using a weak formulation. The weak formulation for maxwell problem is: " find E in H(curl) such as for any E' in H(curl) \int_\Omega -omega/c^2 E.E' + curlE curlE' dx=0 " the first part correspond to the "mass matrix" M, the second to the "stiffness matrix" K you get the mode by taking the eigen vector of (M^-1 K). I have the following code, but unfortunatelly the stiffness matrix doesn't work. I am still trying to find out why. 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)')); M=gf_asm('volumic','M(#1,#1)+=comp(vBase(#1).vBase(#1))(:,k,:,k)', mim,mf); 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); I hope it will help you. Best Regards Louis Kovalevsky On 1 Dec 2013, at 10:18, Jonathan Yik <address@hidden> wrote:
|
[Prev in Thread] | Current Thread | [Next in Thread] |