Re: [deal.II] Vector eigenvalue problem with SLEPc

2021-08-25 Thread Anton Ermakov
>>> *That might not easily be possible with our PETSc interfaces.*





*PETSc has the concept of MatShell, which is in essence what 
deal.II'sLinearOperator is: A class that implements a matrix-vector 
operation,without giving access to individual matrix elements. You could 
trywhether you can implement something derived from MatrixBase that is 
aMatShell operation in the same way as the existing classes derived 
fromMatrixBase model other types of PETSc matrices.*

*I don't know how difficult that would be, though.*

Ok. I see. Seems like using ARPACK might be easier then.

>>> You can do that, but I believe that it's also possible to work with the
original system where one of the two matrices is only semidefinite.

I have redone it with the original system, but it seems that PETSc is 
complaining about missing diagonal entry. 

*[0]PETSC ERROR: Object is in wrong state*

*[0]PETSC ERROR: Matrix is missing diagonal entry 3*

I presume this happens due to a zero block in one (or both) of the original 
matrices. Hmmm, should I then just explicitly see all zero diagonal values 
to 0.0.

Anton.

On Wednesday, August 25, 2021 at 12:20:38 PM UTC-6 Wolfgang Bangerth wrote:

> On 8/25/21 11:53 AM, Anton Ermakov wrote:
> >  I  I¯  A  0  ¯I*eigenvalue +  I¯0 B ¯ I I  * eigenvector =0
> >  I  I_   0 0   _I  I_   C D   _I I
> >  ¯ ¯
> > where eigenvector consists of a vector displacement and a scalar 
> > pressure perturbation:
> > 
> > eigenvector = transpose([U P])
> > 
> > It seems from the literature (e.g., Wang & Bathe 1997) that it is 
> > standard to rewrite this system to eliminate the pressure variable:
> > 
> > eigenvalue * A *  U - B * D^-1 * C * U = 0
> > 
> > or
> > 
> > eigenvalue * A *  U + M * U = 0
> > 
> > with M = - B * D^-1 * C and P = -D^-1 * C * U
>
> You can do that, but I believe that it's also possible to work with the 
> original system where one of the two matrices is only semidefinite.
>
> > What I am having trouble with is how to efficiently represent matrix 
> > products and inverses that go into matrix M. All the matrices (A, B, C, 
> > D) are sparse, but due to the inverse, it seems that M will be a full 
> > matrix and it might be unaffordable to store it in the memory. There was 
> > a discussion on a similar topic here:
> > 
> > https://groups.google.com/g/dealii/c/5P_fv0zg7jg/m/CPYcYrbsDQAJ
> > 
> > It seems that LinearOperator was used in that discussion to manipulate 
> > the matrix products and inverses, which seemed like an approach to 
> > follow if you use ARPACK eigensolvers. However, I am using SLEPc and it 
> > seems the inputs to SLEPc solvers would have to be objects of class 
> > PETScWrappers::MatrixBase 
> > <
> https://nam10.safelinks.protection.outlook.com/?url=https%3A%2F%2Fwww.dealii.org%2Fcurrent%2Fdoxygen%2Fdeal.II%2FclassPETScWrappers_1_1MatrixBase.html&data=04%7C01%7CWolfgang.Bangerth%40colostate.edu%7Ce376b299e56c44690ada08d967f13673%7Cafb58802ff7a4bb1ab21367ff2ecfc8b%7C0%7C0%7C637655108014800275%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C1000&sdata=QQLBoH27O1sZ6GMXDatmucrxAye8edSeN3XNqAybONI%3D&reserved=0>.
>  
>
> > 
> > 
> > So, my question is how to prepare such a matrix M involving matrix 
> > products and matrix inverse efficiently and represent it as 
> > PETScWrappers::MatrixBase 
>
> That might not easily be possible with our PETSc interfaces.
>
> PETSc has the concept of MatShell, which is in essence what deal.II's 
> LinearOperator is: A class that implements a matrix-vector operation, 
> without giving access to individual matrix elements. You could try 
> whether you can implement something derived from MatrixBase that is a 
> MatShell operation in the same way as the existing classes derived from 
> MatrixBase model other types of PETSc matrices.
>
> I don't know how difficult that would be, though.
>
> Best
> W.
>
>
> -- 
> 
> Wolfgang Bangerth email: bang...@colostate.edu
> www: http://www.math.colostate.edu/~bangerth/
>

-- 
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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/88a0aebf-86dd-436f-8f56-8aef2edae4d8n%40googlegroups.com.


Re: [deal.II] Equation implementation

2021-08-25 Thread Toddy Liu
Hi, Daniel

Thank you very much for your valuable advice and I will try with that.

在2021年8月26日星期四 UTC+8 上午5:17:20 写道:

> Toddy,
>
> I would start with thinking about the time discretization scheme you want 
> to use and in particular which terms you want to treat explicitly and which 
> implicitly.
> After you have done that, you might still end up with a nonlinear problem 
> in which case you should have a look at 
> https://www.dealii.org/current/doxygen/deal.II/step_15.html and 
> https://www.dealii.org/current/doxygen/deal.II/step_25.html.
>
> Best,
> Daniel
>
> Am Mi., 25. Aug. 2021 um 09:11 Uhr schrieb Toddy Liu  >:
>
>> Dear Deal.II community,
>>
>> I'm getting trouble in constructing system matrix of a complex PDE. My 
>> PDE is a time-dependent equation, and in simplicity it seems like 
>>
>> [image: screenshot.png]
>> U is the variable, t denotes the time. D, alpha, beta are constant. So 
>> could you please give me some advice about doing this or any reference in 
>> existing dealii examples.
>>
>> Thank you very much.
>>
>> -- 
>> 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 dealii+un...@googlegroups.com.
>> To view this discussion on the web visit 
>> https://groups.google.com/d/msgid/dealii/b239ae4a-a114-47a7-ad59-a275f8771239n%40googlegroups.com
>>  
>> 
>> .
>>
>

-- 
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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/77c53750-7e17-402b-a927-71891a8388a1n%40googlegroups.com.


Re: [deal.II] Equation implementation

2021-08-25 Thread Daniel Arndt
Toddy,

I would start with thinking about the time discretization scheme you want
to use and in particular which terms you want to treat explicitly and which
implicitly.
After you have done that, you might still end up with a nonlinear problem
in which case you should have a look at
https://www.dealii.org/current/doxygen/deal.II/step_15.html and
https://www.dealii.org/current/doxygen/deal.II/step_25.html.

Best,
Daniel

Am Mi., 25. Aug. 2021 um 09:11 Uhr schrieb Toddy Liu :

> Dear Deal.II community,
>
> I'm getting trouble in constructing system matrix of a complex PDE. My PDE
> is a time-dependent equation, and in simplicity it seems like
>
> [image: screenshot.png]
> U is the variable, t denotes the time. D, alpha, beta are constant. So
> could you please give me some advice about doing this or any reference in
> existing dealii examples.
>
> Thank you very much.
>
> --
> 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 dealii+unsubscr...@googlegroups.com.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/dealii/b239ae4a-a114-47a7-ad59-a275f8771239n%40googlegroups.com
> 
> .
>

-- 
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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/CAOYDWb%2Bd7hp%2Bww7w9AKmYm1JnSwt8ZAFmig%3D5W_%3DZ-O6R14ORg%40mail.gmail.com.


Re: [deal.II] Vector eigenvalue problem with SLEPc

2021-08-25 Thread Wolfgang Bangerth

On 8/25/21 11:53 AM, Anton Ermakov wrote:

  I  I¯  A  0  ¯I    *eigenvalue     +  I¯    0 B ¯ I I  * eigenvector =0
  I  I_   0 0   _I                                  I_   C D   _I I
  ¯                                                                     ¯
where eigenvector consists of a vector displacement and a scalar 
pressure perturbation:


eigenvector = transpose([U P])

It seems from the literature (e.g., Wang & Bathe 1997) that it is 
standard to rewrite this system to eliminate the pressure variable:


eigenvalue * A *  U - B * D^-1 * C * U = 0

or

eigenvalue * A *  U + M * U = 0

with M = - B * D^-1 * C and P = -D^-1 * C * U


You can do that, but I believe that it's also possible to work with the 
original system where one of the two matrices is only semidefinite.


What I am having trouble with is how to efficiently represent matrix 
products and inverses that go into matrix M. All the matrices (A, B, C, 
D) are sparse, but due to the inverse, it seems that M will be a full 
matrix and it might be unaffordable to store it in the memory. There was 
a discussion on a similar topic here:


https://groups.google.com/g/dealii/c/5P_fv0zg7jg/m/CPYcYrbsDQAJ

It seems that LinearOperator was used in that discussion to manipulate 
the matrix products and inverses, which seemed like an approach to 
follow if you use ARPACK eigensolvers. However, I am using SLEPc and it 
seems the inputs to SLEPc solvers would have to be objects of class 
PETScWrappers::MatrixBase 
. 



So, my question is how to prepare such a matrix M involving matrix 
products and matrix inverse efficiently and represent it as 
PETScWrappers::MatrixBase 


That might not easily be possible with our PETSc interfaces.

PETSc has the concept of MatShell, which is in essence what deal.II's 
LinearOperator is: A class that implements a matrix-vector operation, 
without giving access to individual matrix elements. You could try 
whether you can implement something derived from MatrixBase that is a 
MatShell operation in the same way as the existing classes derived from 
MatrixBase model other types of PETSc matrices.


I don't know how difficult that would be, though.

Best
 W.


--

Wolfgang Bangerth  email: bange...@colostate.edu
   www: http://www.math.colostate.edu/~bangerth/

--
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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/9753bb31-5c06-88c5-036f-008196d1d472%40colostate.edu.


[deal.II] Vector eigenvalue problem with SLEPc

2021-08-25 Thread Anton Ermakov
Dear Deal.II community,

I am working on implementing a vector eigenproblem in deal.ii, starting 
from the step-36 example. The question I am studying is the normal 
oscillation of a fluid planet. The eigenproblem is of this type:
 _  _
 I  I¯  A  0  ¯I*eigenvalue +  I¯0 B ¯ I I  * eigenvector =0
 I  I_   0 0   _I  I_   C D   _I I
 ¯ ¯
where eigenvector consists of a vector displacement and a scalar pressure 
perturbation: 

eigenvector = transpose([U P])

It seems from the literature (e.g., Wang & Bathe 1997) that it is standard 
to rewrite this system to eliminate the pressure variable:

eigenvalue * A *  U - B * D^-1 * C * U = 0

or 

eigenvalue * A *  U + M * U = 0

with M = - B * D^-1 * C and P = -D^-1 * C * U 

What I am having trouble with is how to efficiently represent matrix 
products and inverses that go into matrix M. All the matrices (A, B, C, D) 
are sparse, but due to the inverse, it seems that M will be a full matrix 
and it might be unaffordable to store it in the memory. There was a 
discussion on a similar topic here:

https://groups.google.com/g/dealii/c/5P_fv0zg7jg/m/CPYcYrbsDQAJ

It seems that LinearOperator was used in that discussion to manipulate the 
matrix products and inverses, which seemed like an approach to follow if 
you use ARPACK eigensolvers. However, I am using SLEPc and it seems the 
inputs to SLEPc solvers would have to be objects of class 
PETScWrappers::MatrixBase 

. 

So, my question is how to prepare such a matrix M involving matrix products 
and matrix inverse efficiently and represent it as PETScWrappers::MatrixBase 

 suitable 
for a SLEPc eigensolver. 

In addition, I am interested in the interior eigenvalues. It seems that 
shift-invert spectral transformation is available in SLEPc - this is what I 
plan to use. However, SLEPc also has a polynomial filtering spectral 
transformation. It does not seem that deal.ii provides a wrapper to that 
one though. So, I wonder if anyone has had expertise with making an 
interface in deal.ii to use SLEPc's polynomial filtering transformation.

Thanks a lot,
Anton.



-- 
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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/426e3d15-a42f-4e90-a3f7-117c66e71c92n%40googlegroups.com.


[deal.II] Equation implementation

2021-08-25 Thread Toddy Liu
Dear Deal.II community,

I'm getting trouble in constructing system matrix of a complex PDE. My PDE 
is a time-dependent equation, and in simplicity it seems like 

[image: screenshot.png]
U is the variable, t denotes the time. D, alpha, beta are constant. So 
could you please give me some advice about doing this or any reference in 
existing dealii examples.

Thank you very much.

-- 
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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/b239ae4a-a114-47a7-ad59-a275f8771239n%40googlegroups.com.