Dear Johnathan,
The expressions I gave to you are correct, you can enforce the symmetry
> M=gf_asm('volumic','M(#1,#1)+=sym(comp(vBase(#1).vBase(#1))(:,k,:,k))',
> mim,mf);
>
> K=gf_asm('volumic','M(#1,#1)+=sym(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);
Then you need to calculate the eigen vectors of M^-1 K
Be aware that the Kernel is very big.
I cite Ronan:
"Your matrix K should be positive semi-definite with a large kernel equal to
the number of nodes in your mesh minus one (modulo the boundary conditions)."
Best Regards
Louis
On 1 Dec 2013, at 15:43, Louis Kovalevsky <[email protected]> wrote:
> 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 <[email protected]> wrote:
>
>> Hello,
>>
>> I have been trying to calculate the eigenmodes of light travelling through a
>> dielectric waveguide. To do this, I have been using the assembly tools for
>> arbitrary matrices. However, when I attempt to solve for eigenvectors of
>> the stiffness matrix using ARPACK, I do not get sensible results.
>>
>> I have been attempting to generate matrices for the Maxwell Wave equation
>> using:
>>
>> Del^2 Psi - (omega^2 / c^2 ) Psi
>>
>> For the first term, I simply use the laplacian given in the getfem
>> documentation: M(#1, #1)+=comp(Grad(#1).Grad(#1))(:,k,:,k). However, I am
>> unsure as to how I should generate the second term. However, I am not sure
>> how to generate the second term in the equation. Previously, I have used:
>> (comp(Base(#1).Base(#1).Base(#2))(:,:,j).h(j)). Another technique I have
>> considered was to take the vector generated by interpolating the data set
>> representing (omega^2 / c^2 ), converting the elements of the vector to
>> elements of a diagonal matrix, and simply adding the result to the matrix of
>> the first term. Could anyone tell me if either of these approaches is
>> correct, or if the correct procedure is some method that I have not yet
>> thought of?
>>
>> Yours,
>>
>> Johnathan Yik
>> _______________________________________________
>> 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