open TELEMAC-MASCARET The mathematically superior suite of solvers

Hello everyone,

I am trying to solve non-linear partial differential equation coupled with Telemac 2d. I have two builds two lumped mass matrices as shown in the attached picture.

I would appreciate if you can help in solving this,

Note,

If I build a mass matrix as shown at the bottom of the image, then multiply by vector then call lump routine, does represent the same as I posted at the top of the image matrix?

How about the second Matrix, Dij, I have gradPSIj.PSIJ .gradPSIj, how to build this one?

Thanks,
AMANJ

Plz confirm this too!

Hello,

This last is confirmed only if you add a sum on j on the left hand side.

In the first post it seems to me that terms S and C, being nodal values and not functions, can get out of the integral. Then with the Kronecker symbol you are only left with the diagonal of the mass-matrix multiplied by a vector, this will do Oij.

The other matrix is more difficult, it is close to the matrices done for SUPG. However my general impression is that something is fishy somewhere, for example Oij does not looks like the result of a variational formulation. Maybe showing us the original equations that you want to solve, in the continuum would be helpful.

To be more specific, when you multiply your matrix by a vector, and sum all the components of the vector, using the fact that the summ of all basis makes 1 everywhere, it must give you the integral of something physical, like a mass or a momentum, not something that still depends on the basis.

With best regards,

Jean-Michel Hervouet

Hello dear Jean-Michel Hervouet,

Thank you for always support and help,

Here my original equation and weak statement.



S, C and Kr are nodal values and each one is a vector with size 1:NPOIN and each one need to interpolate using shape function to find the element value.

I need your help in order to use BIEF Library in order to build matrices.

In my case, for building O, which is lumped-mass matrix, need first to build MATMAS
SUCH AS THIS ONE IN TELEMAC2D:
FORMUL=MATMAS
CALL MATRIX(AM1,'M=N ',FORMUL,IELMH,IELMH,
& SL1,TE3,S,S,S,S,S,MESH,MSK,MASKEL)

THEN I HAVE TO LUMPED SUCH AS THIS EXAMPLE:

CALL LUMP(T2,AM1,MESH,AGGLOU)

THEN MUTIPLY BY VECTOR WHAT I HAVE SUCH AS :

(sS0+NU*C) WHERE s AND CARE NODAL VALUE AND S0 AND NU ARE CONSTANT.

DOES THIS WAY SOLVE MY Oij OR NOT?

FOR BUILDING Dij, NEED TO MODIFIY ONE OF THE INTEGRATION MATRICES OR CAN I USE MEAN INTERPOLATION FOR EXAMPLE ( KR ELEMENT=(KR1+KR2+KR2)/3 AND MULTIPLY THE MATRIX Dij, SAME FOR s AND C INSTEAD OF USING INTERPOLATION FUNCTION.

THANKS,

AMANJ

Hello,

If you need a fully lumped matrix, then it is a vector, so rather than really building it and then lumping it, you should rather do that on a paper, i.e. summing O(i,j) over j, to get the vector Olumped(i). Generally terms disappear and it simplifies using the fact that the sum of all test functions is equal to 1. For example lumping a diffusion matrix gives zero. Any function multiplied by gradient(PSIj) would cancel also. In this respect Kr in your basic equation is strange because it seems to be already a discretisation, as you say it is a nodal value.
You can always simplify, like taking an average Kr on an element, but beware, it must not spoil mass conservation or other important properties. What we do for example in friction terms is that the friction coefficient within an integral (as a continuous function) is then put outside the integral, choosing a nodal value. There is a theorem saying that there exist an average value of the friction for which we have an equality. We do not know the average value but taking the nodal value is not bad, especially when the friction is locally constant...

With best regards,

Jean-Michel Hervouet

I have questions related to 3D linear triangular prism,
first of all, the period of integration: is it same as the following for simple PSI dPSI or mass matrix:

int(int(int((1-ALPHA-BETA)*(1-GAMA),0..1-BETA),0..1),-1..1) or in PSI1*PSI1*dPSI
int(int(int((1-ALPHA-BETA)^2*(1-GAMA)^,0..1-BETA),0..1),-1..1) or something different?

second about JACOBIAN matrix for prism, its still function of ALPHA,BETA, and GAMA,

did you use centroid interpolation to vanish them or something else, someone used 1/3 for ALPHA and BETA and 0 for GAMA and someone else 2/3 instead of 1/3.

I would like to know how the integration performed.

I tried Maple for integration but I didnt get same coefficents what we have for example in file MT00PP.

Thanks,

AMANJ

Hello,

I see that you are looking at the core. The period of integration seems OK and yes , you are right, there is a problem with diffusion-like matrices, there remains a Jacobian in the denominator. This topic is treated page 211 and 212 of my book, I just summarise here :

* you can take an average Jacobian in the denominator.

* In the case of the diffusion matrix, it must be in the end a M-matrix, that is to say that the diagonal has only negative values and the off-diagonal terms positive values (this ensures monotony). If you just compute the integrals following the theory of finite elements it will not give a M-matrix. The simplification that you can do is to compute the vertical derivatives using only points on the same vertical.

* always keep the good properties of the true matrix: like sum of lines =0 and/or sum of columns=0. For the diffusion matrix, we compute the extra-diagonal terms and then the diagonal is deduced by the fact that the summ of columns is 0.

I hope this helps, I remember spending a lot of time on this,

With best regards,

Jean-Michel Hervouet

Hello dear Jean-Michel Hervouet ,

I tried to understand your explanation but unfortunately I couldn't. Is there any reference explaining this in detail by showing an example how to calculate the coefficients of diffusion matrix?

I have my jacobian |Je|=h1+(h2-h1)*alpha+(h3-h1)*Beta,

In the matrix 3x3 and vector 3x1 I have alpha and beta and as you said its impossible analytically to integrate.

I didn't get what did you mean or how to take average of Jaccobian?

which lines or coloumns you mean equal to 0?

Is it possible if you share details of calculating diffusion matrix with us?

I am really struggling with this.

Thanks,

AMANJ

Hello,

I have little more than what is explained in my book and the implementation in sources. Actually to take an average Jacobian you can just choose the value at a point in the middle of the prism (it makes a big simplification). The rows and columns are the rows and columns of the matrix. If you have a gradient(PSI(i)) in a variational formulation with nothing depending on a specific point in the rest of the formula, the sum over all points i, i.e. the sum of a column of the matrix, will be 0 because the sum of all PSI(I) is one if PSI(i) are the test functions at points i. A diffusion matrix is grad(PSI(i))*grad(PSI(j)) so the sum of rows and columns is 0. If you take a diffusion matrix on triangles, this is on every element a priori 9 coefficients to compute, with symmetry it is actually only 6, using the "magic square" properties : sum of rows and columns equal to 0, it is in fact only 3 coefficients to compute. You can first use these properties to check your computation (of course if your matrix has such properties).

With best regards,

Jean-Michel Hervouet