Discussion Forum

Hi Marc,

I'm working on a plate vibration problem using eigenvalue analysis. The dependent variables are u, v, w and arx, ary, arz. I export the mode profiles sampled at a grid and for some reason the modes are not orthogonal.

I would first recommend not exporting the values but perform the orthogonality check directly in COMSOL, which reduces the sources of possible errors and will be more accurate because COMSOL will integrate with knowledge of interpolation functions and hence more accurate than based on samples.
I encourage you to send us the question and a suitable model in regular support

support@comsol.com

With best regards from Technopark Zurich.

Sven Friedel

Thank you! I will send an e-mail to the support you provided me.

As a side note, If I integrate the eigenfunctions in Comsol using with(i,u)*with(j, u), ... I get very similar results, but it's a very good suggestion to remove the integration error.

The displacements uvw and the rotations ar have different physical dimensions. So summing the expression like you suggest will not give anything useful.

Let 'U' be the physical displacement, and 'u', 'a' be the degrees of freedom. Then the mathematical inner product of two eigenmodes is



where the mass density has been assumed to be constant.

Functions which are odd in integrates to zero over the thickness 'd', so we end up with



Now, for a normal plate bending problem, only the transverse displacement is non-zero. For a thin plate the rotations can be ignored, but they will give a correction term for a thick plate. The correction is however of the order of



where L is a characteristic width of the plate. So the plate must be rather thick to make it worthwhile to incorporate the rotations. But if you are checking orthogonality, you may still come several orders of magnitudes closer to zero by incorporating the rotations.

If you are going to perform these integrations, it is a good idea to use 'Scaling of eigenvectors' = 'Mass matrix', because then you can check that you get exactly 1.0 when you integrate a mode against itself.

Regards,
Henrik

Thank you very much for your comment. That was a big mistake in our part.

I modified our model so it calculates the surface integral of:

with(21,u)*with(22,u)+with(21,v)*with(22,v)+with(21,w)*with(22,w)+(with(21,arx)*with(22,arx)+with(21,ary)*with(22,ary)+with(21,arz)*with(22,arz))*((plate.d)^2)/12

And selected the mass matrix normalization.

For each mode I also calculate:

u*u+v*v+w*w+(arx*arx+ary*ary+arz*arz)*((plate.d)^2)/12

However the results are still problematic. The product of the modes by themselves is 0.63694 instead of 1, and the product of one mode by another is 0.0518 (~8% "non-orthogonality"). Most of the latter comes from the term w*w (removing the rest affects only the last significant digit.

Any suggestion on what may be going on?

Hi,

You are missing a factor plate.rho*plate.d in front of the whole expression (see the integrals in my previous posting).

The bad orthogonality is trickier to guess though. But it is a high mode number you are working with, and higher modes have lower accuracy. Try

a) Inner product between two of the lower modes to get a comparison
b) Finer mesh

Regards,
Henrik

Thank you, this is very helpful. Now the modes integrate to 1 except two of them which are a pair of complex-conjugate values. I expect the matrix to be symmetric, so I'm very surprised about that. The complex modes are the 27 and 26, which I don't think are "that high" for a plate. Keep in mind that most modes have a 2-fold degeneracy due to symmetry.

The modes I'm interested in are the 21 and 22. The first 6 modes are rigid body motions. I'm using a spring foundation of 100 N/m total stiffness to stabilize the plate. So I'm looking at the modes 15-16. The mesh has 40k elements.

I attach my model for reference purposes.

Marc

Hi Marc,

This is also a support case, but since we already have this forum discussion, I will answer it here. Part of the answers have a general interest.

1. With that many elements, your mesh density is more than sufficient for what you would need, so that is not a problem.

2. Complex results can sometimes be seen in eigenvalue problems where you, for physics reasons, would not expect it. The eigenvalue solver treats any eigenvalue problem as complex, so retrieving a real valued solution relies on the numerical precision of the problem. There is no assumption made about real valued eigenfrequencies and eigenmodes.

3. The lack of orthogonality between these two modes is related to lack of numerical precision This is sometimes a problem for pairs of degenerate eigenmodes. For other pairs there is more or less perfect orthogonality.

In your case, there are three things that are sensitive from the numerical point of view, and which can affect both item 2) and 3):

a) There are rigid body modes. Actually, the springs you add are not necessary; COMSOL can handle unconstrained eigenvalue problems. But you should set 'Search for eigenfrequencies around' in the eigenfrequency solver. See the quotes from the user's guide below:

"It is possible to compute eigenfrequencies for structures which are not fully constrained; this is sometimes referred to as free-free modes. For each possible rigid body mode, there is one eigenvalue which in theory is zero."

"In practice, the natural frequencies are not computed as exactly zero, but can appear as small numbers which can even be negative or complex. If rigid body modes are present in the model, then it is important to use a nonzero value in the Search for eigenfrequencies around text field in the settings for the Eigenfrequency study step. The value should reflect the order of magnitude of the first important nonzero eigenfrequency."

b) There are degenerate eigenfrequencies. The orthogonalization is much more sensitive then. I notice that you have already used a tight tolerance for the eigenvalue solver. Some suggestion for improving the accuracy further:

- You can tighten the tolerance even more. Check in the solver log that your changed setting acturally forces one more iteration. If it looks like below
Iter ErrEst Nconv
1 2e-031 30
then you can set tolerance to 1e-32.

- You can increase the accuracy of the eigenvalue orthogonalization by setting "Dimension of Krylov space" to something high, say three times the number of eigenmodes you are computing.

- The most efficient remedy is however to disturb the mesh slightly. If you mesh the square with 199x200 elements instead of 200x200 elements, the orthogonality between the modes is restored. The reason is that the perfect duplicity of the eigenvalues is removed. The eigenvalues will still be reported as equal, though.

c) You are carrying along a lot of unnecessary degrees of freedom. If it is a pure plate bending problem, you can decrease the size of the problem significantly by deselecting 'Use 3D formulation' at the settings for the physics interface. Note that you will then only have three rigid-body modes instead of six, so the inner product should be (18,19) instead of (21,22).

Regards,
Henrik