Dear Wolfgang:
Thank you for your answer. I have checked the matrices and I printed the
case for a 2 by 2 mesh from both Matlab and deal.II. With a simple change
of rows and columns, they are exactly the same ( see the files attached).
Not sure yet why they result in different results of the eigenvalue
problem. There is a very importent detail that I want to point out here.
When I solved on the very course mesh (2 by 2) I got the following error,
which is not the case when I solved on a 4 by 4 mesh for example.
An error occurred in line <791> of file
</usr/local/include/deal.II/lac/arpack_solver.h> in function
void dealii::ArpackSolver::solve(const MatrixType1&, const
MatrixType2&, const INVERSE&, std::vector<std::complex<double> >&,
std::vector<ValueType>&, unsigned int) [with VectorType =
dealii::Vector<double>; MatrixType1 =
dealii::LinearOperator<dealii::Vector<double>, dealii::Vector<double>,
dealii::internal::LinearOperatorImplementation::EmptyPayload>; MatrixType2
= dealii::LinearOperator<dealii::Vector<double>, dealii::Vector<double>,
dealii::internal::LinearOperatorImplementation::EmptyPayload>; INVERSE =
dealii::SparseDirectUMFPACK]
The violated condition was:
false
Additional information:
Error with dnaupd, info -2. Check documentation of ARPACK
Dear Mathias:
Thank you for your answer. The matrices are exactly the same in both codes.
The different is that I solve for the recipricals in deal.II. To be precise,
let AA= Bt A^inv B
and let M be a mass matrix
where B, A, M are the same matrices in both Matlab and deal.II.
In Matlab: I solve eig(AA, M) (AA x = \lambda M x) then I seeks the
square root of the smalles \lambda
In dealii: I solve using ArpackSolver solving M x = 1/\Lambda AA x)
so I do first
M_inv_umfpack.initialize(M);
auto M = linear_operator(M);
auto M_inv = linear_operator(M, M_inv_umfpack);
then, I create a linear operator of the other matrices and the inverse of A
and perform AA= Bt * A^inv * B
lastly I solve : arpack_solver.solve( M, AA, M_inv_umfpack,
eigenvalues,eigenvectors);
So I am seeking the largest eigenvalue and \lambda = 1/largest eigenvalue
and the infsup constant is the square root of that.
Could solving for the recipricals is the issue? if yes what should I fill
in the third intery of the solving if Iam doing
arpack_solver.solve( AA, M, ????????, eigenvalues,eigenvectors);
One more time, your questions and help are highly appreciated.
Best,
Najwa
On Friday, May 31, 2024 at 12:40:04 AM UTC+3 Matthias Maier wrote:
> Do you use the same "mass matrix" in both cases, meaning do you call
> eig(A, M) where A is the stiffness matrix and M is the mass matrix in
> matlab?
>
>
> You will need to compute the spectrum of the eigenvalue problem
>
> A x = \lambda M x
>
>
> Best,
> Matthias
>
>
> On Thu, May 30, 2024, at 10:55 CDT, Najwa Alshehri <[email protected]>
> wrote:
>
> > Dear all,
> >
> >
> > Thank you for your always support and help.
> >
> >
> > I am encountering an issue while solving an eigenvalue problem related
> to
> > the computation of the discrete inf-sup constant for a saddle point
> > problem. Specifically, I am solving the following system:
> >
> >
> > Bt A^-1 B eigenvector = eigenvalue M eigenvector,
> >
> >
> > where Bt A^-1 B is clearly symmetric and at least semi-positive
> definite,
> > and M is a symmetric positive definite matrix.
> >
> >
> > The problem arises when comparing the inf-sup constants obtained using
> > Matlab and dealii (solving on the same mesh ):
> >
> >
> > 1. Matlab: Using the "eig" solver, the inf-sup constant decays to zero
> > as the mesh is refined.
> > 2. deal.II: Using the "ArpackSolver", the inf-sup constant is bounded
> > away from zero.
> >
> > The discrepancy is puzzling, as both methods should give consistent
> > results. To illustrate this issue, I have attached a figure comparing
> the
> > two results and included in this post the two data files obtained from
> > Matlab and deal.II.
> >
> >
> > I would greatly appreciate any insights or suggestions on what might be
> > causing the discrepancy.
> >
> >
> > Your help is highly appreciated,
> >
> > Najwa
>
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/ed61d756-3f6c-4d1c-a7ae-7337def7ebc6n%40googlegroups.com.
Mu2 =
0.7229 -0.1385 0 -0.1385 -0.3193 0 0
0 0
-0.1385 1.4459 -0.1385 -0.3193 -0.2771 -0.3193 0
0 0
0 -0.1385 0.7229 0 -0.3193 -0.1385 0
0 0
-0.1385 -0.3193 0 1.4459 -0.2771 0 -0.1385
-0.3193 0
-0.3193 -0.2771 -0.3193 -0.2771 2.8918 -0.2771 -0.3193
-0.2771 -0.3193
0 -0.3193 -0.1385 0 -0.2771 1.4459 0
-0.3193 -0.1385
0 0 0 -0.1385 -0.3193 0 0.7229
-0.1385 0
0 0 0 -0.3193 -0.2771 -0.3193 -0.1385
1.4459 -0.1385
0 0 0 0 -0.3193 -0.1385 0
-0.1385 0.7229
C2 =
0.1266 0.1266 0 0.1266 0.1266 0 0
0 0
0 0.1266 0.1266 0 0.1266 0.1266 0
0 0
0 0 0 0.1266 0.1266 0 0.1266
0.1266 0
0 0 0 0 0.1266 0.1266 0
0.1266 0.1266
Nlambda =
0.5130 0 0 0
0 0.5130 0 0
0 0 0.5130 0
0 0 0 0.5130
rearranging:
C2 =
0.1266 0.1266 0.1266 0.1266 0 0 0
0 0
0 0.1266 0 0.1266 0.1266 0.1266 0
0 0
0 0 0.1266 0.1266 0 0 0.1266
0.1266 0
0 0 0 0.1266 0 0.1266 0
0.1266 0.1266
Mu2 =
0.7229 -0.1385 -0.1385 -0.3193 0 0 0
0 0
-0.1385 1.4459 -0.3193 -0.2771 -0.1385 -0.3193 0
0 0
-0.1385 -0.3193 1.4459 -0.2771 0 0 -0.1385
-0.3193 0
-0.3193 -0.2771 -0.2771 2.8918 -0.3193 -0.2771 -0.3193
-0.2771 -0.3193
0 -0.1385 0 -0.3193 0.7229 -0.1385 0
0 0
0 -0.3193 0 -0.2771 -0.1385 1.4459 0
-0.3193 -0.1385
0 0 -0.1385 -0.3193 0 0 0.7229
-0.1385 0
0 0 -0.3193 -0.2771 0 -0.3193 -0.1385
1.4459 -0.1385
0 0 0 -0.3193 0 -0.1385 0
-0.1385 0.7229M_u2 :
0.7229 -0.1385 -0.1385 -0.3193 0.0000 0.0000 0.0000
0.0000 0.0000
-0.1385 1.4459 -0.3193 -0.2771 -0.1385 -0.3193 0.0000
0.0000 0.0000
-0.1385 -0.3193 1.4459 -0.2771 0.0000 0.0000 -0.1385
-0.3193 0.0000
-0.3193 -0.2771 -0.2771 2.8918 -0.3193 -0.2771 -0.3193
-0.2771 -0.3193
0.0000 -0.1385 0.0000 -0.3193 0.7229 -0.1385 0.0000
0.0000 0.0000
0.0000 -0.3193 0.0000 -0.2771 -0.1385 1.4459 0.0000
-0.3193 -0.1385
0.0000 0.0000 -0.1385 -0.3193 0.0000 0.0000 0.7229
-0.1385 0.0000
0.0000 0.0000 -0.3193 -0.2771 0.0000 -0.3193 -0.1385
1.4459 -0.1385
0.0000 0.0000 0.0000 -0.3193 0.0000 -0.1385 0.0000
-0.1385 0.7229
C2 :
0.1266 0.1266 0.1266 0.1266 0.0000 0.0000 0.0000
0.0000 0.0000
0.0000 0.1266 0.0000 0.1266 0.1266 0.1266 0.0000
0.0000 0.0000
0.0000 0.0000 0.1266 0.1266 0.0000 0.0000 0.1266
0.1266 0.0000
0.0000 0.0000 0.0000 0.1266 0.0000 0.1266 0.0000
0.1266 0.1266
M_lambda :
0.5130 0.0000 0.0000 0.0000
0.0000 0.5130 0.0000 0.0000
0.0000 0.0000 0.5130 0.0000
0.0000 0.0000 0.0000 0.5130
