Re: [petsc-users] question about MatSetLocalToGlobalMapping

[Zhang, Hong via petsc-users]

Pierre,


1) Is there any hope to get PDIPDM to use a MatNest?

KKT matrix is indefinite and ill-conditioned, which must be solved using a 
direct matrix factorization method.
But you are using PCBJACOBI in the paper you attached?
In any case, there are many such systems, e.g., a discretization of Stokes 
equations, that can be solved with something else than a direct factorization.

The KKT in this paper is very special, each block represents a time step with 
loose and weak couplings. We tested larger and slightly varying physical 
parameters, PCZBJACOBI fails convergence while CHOLESKY encounters out of 
memory.
For the current implementation, we use MUMPS Cholesky as default. To use 
MatNest, what direct solver to use, SCHUR_FACTOR? I do not know how to get it 
work.
On the one hand, MatNest can efficiently convert to AIJ or SBAIJ if you want to 
stick to PCLU or PCCHOLESKY.
On the other hand, it allows to easily switch to PCFIELDSPLIT which can be used 
to solve saddle-point problems.

MatNest would be a natural way for these type of applications. PETSc has 
MatConvert_Nest_AIJ(), which could enable users to assemble their matrices in 
MatNest, and apply mumps/superlu direct solvers. I'm not concerned about the 
overhead of matrix conversion, because matrix factorization dominates 
computation in general.
2) Is this fixed 
https://lists.mcs.anl.gov/pipermail/petsc-dev/2020-September/026398.html ?
I cannot get users to transition away from Ipopt because of these two missing 
features.

The existing pdipm is the result of a MS student intern project. None of us 
involved are experts on the optimization solvers. We made a straightforward 
parallelization of Ipopt. It indeed needs further work, e.g., more features, 
better matrix storage, convergence criteria... To our knowledge, parallel pdipm 
is not available other than our pdipm.
Ipopt can use MUMPS and PARDISO internally, so it’s in some sense parallel 
(using shared memory).
Also, this is not a very potent selling point.
My users that are satisfied with Ipopt as a "non-parallel" black box don’t want 
to have to touch part of their code just to stick it in a parallel black box 
which is limited to the same kind of linear solver and which has severe 
limitations with respect to Hessian/Jacobian/constraint distributions.

We saw exiting 'parallel' pdipm papers, which conduct sequential computation 
and send KKT matrices to multiple processors, then they show scaling results on 
matrix computation only. Again, I wish someone can help to improve tao/pdipm. 
This is a useful solver.
We should improve our pdipm.
Hong

On 20 Apr 2023, at 5:47 PM, Zhang, Hong via petsc-users 
mailto:petsc-users@mcs.anl.gov>> wrote:

Karthik,
We built a KKT matrix in  TaoSNESJacobian_PDIPM() (see 
petsc/src/tao/constrained/impls/ipm/pdipm.c) which assembles several small 
matrices into a large KKT matrix in mpiaij format. You could take the same 
approach to insert P and P^T into your K.
FYI, I attached our paper.
Hong


From: petsc-users 
mailto:petsc-users-boun...@mcs.anl.gov>> on 
behalf of Matthew Knepley mailto:knep...@gmail.com>>
Sent: Thursday, April 20, 2023 5:37 AM
To: Karthikeyan Chockalingam - STFC UKRI 
mailto:karthikeyan.chockalin...@stfc.ac.uk>>
Cc: petsc-users@mcs.anl.gov<mailto:petsc-users@mcs.anl.gov> 
mailto:petsc-users@mcs.anl.gov>>
Subject: Re: [petsc-users] question about MatSetLocalToGlobalMapping

On Thu, Apr 20, 2023 at 6:13 AM Karthikeyan Chockalingam - STFC UKRI via 
petsc-users mailto:petsc-users@mcs.anl.gov>> wrote:
Hello,

I created a new thread, thought would it be more appropriate (and is a 
continuation of my previous post). I want to construct the below K matrix 
(which is composed of submatrices)

K = [A P^T
   P   0]

Where K is of type MatMPIAIJ. I first constructed the top left [A] using 
MatSetValues().

Now, I would like to construct the bottom left [p] and top right [p^T] using 
MatSetValuesLocal().

To use  MatSetValuesLocal(),  I first have to create a local-to-global mapping 
using ISLocalToGlobalMappingCreate. I have created two mapping row_mapping and 
column_mapping.

I do not understand why they are not the same map. Maybe I was unclear before. 
It looks like you have two fields, say phi and lambda, where lambda is a 
Lagrange multiplier imposing some constraint. Then you get a saddle point like 
this. You can imagine matrices

  (phi, phi)--> A
  (phi, lambda) --> P^T
  (lambda, phi) --> P

So you make a L2G map for the phi field and the lambda field. Oh, you are 
calling them row and col map, but they are my phi and lambda
maps. I do not like the row and col names since in P they reverse.

Q1) At what point should I declare MatSetLocalToGlobalMapping – is it just 
before I use MatSetValuesLocal()?

Okay, it is good you are asking this because my thinking was somewhat confused. 
I think the pr

Re: [petsc-users] question about MatSetLocalToGlobalMapping

[Karthikeyan Chockalingam - STFC UKRI via petsc-users]

Thank you Hong. I will look into it.

Best,
Karthik.

From: Zhang, Hong 
Date: Thursday, 20 April 2023 at 16:47
To: Matthew Knepley , Chockalingam, Karthikeyan (STFC,DL,HC) 

Cc: petsc-users@mcs.anl.gov 
Subject: Re: [petsc-users] question about MatSetLocalToGlobalMapping
Karthik,
We built a KKT matrix in  TaoSNESJacobian_PDIPM() (see 
petsc/src/tao/constrained/impls/ipm/pdipm.c) which assembles several small 
matrices into a large KKT matrix in mpiaij format. You could take the same 
approach to insert P and P^T into your K.
FYI, I attached our paper.
Hong

From: petsc-users  on behalf of Matthew 
Knepley 
Sent: Thursday, April 20, 2023 5:37 AM
To: Karthikeyan Chockalingam - STFC UKRI 
Cc: petsc-users@mcs.anl.gov 
Subject: Re: [petsc-users] question about MatSetLocalToGlobalMapping

On Thu, Apr 20, 2023 at 6:13 AM Karthikeyan Chockalingam - STFC UKRI via 
petsc-users mailto:petsc-users@mcs.anl.gov>> wrote:

Hello,



I created a new thread, thought would it be more appropriate (and is a 
continuation of my previous post). I want to construct the below K matrix 
(which is composed of submatrices)



K = [A P^T

   P   0]



Where K is of type MatMPIAIJ. I first constructed the top left [A] using 
MatSetValues().



Now, I would like to construct the bottom left [p] and top right [p^T] using 
MatSetValuesLocal().



To use  MatSetValuesLocal(),  I first have to create a local-to-global mapping 
using ISLocalToGlobalMappingCreate. I have created two mapping row_mapping and 
column_mapping.

I do not understand why they are not the same map. Maybe I was unclear before. 
It looks like you have two fields, say phi and lambda, where lambda is a 
Lagrange multiplier imposing some constraint. Then you get a saddle point like 
this. You can imagine matrices

  (phi, phi)--> A
  (phi, lambda) --> P^T
  (lambda, phi) --> P

So you make a L2G map for the phi field and the lambda field. Oh, you are 
calling them row and col map, but they are my phi and lambda
maps. I do not like the row and col names since in P they reverse.


Q1) At what point should I declare MatSetLocalToGlobalMapping – is it just 
before I use MatSetValuesLocal()?

Okay, it is good you are asking this because my thinking was somewhat confused. 
I think the precise steps are:

  1) Create the large saddle point matrix K

1a) We must call 
https://petsc.org/main/manualpages/Mat/MatSetLocalToGlobalMapping/ on it. In 
the simplest case, this just maps
   the local rows numbers [0, Nrows) to the global rows numbers 
[rowStart, rowStart + Nrows).

  2) To form each piece:

2a) Extract that block using 
https://petsc.org/main/manualpages/Mat/MatGetLocalSubMatrix/

   This gives back a Mat object that you subsequently restore using 
https://petsc.org/main/manualpages/Mat/MatRestoreLocalSubMatrix/

 2b) Insert values using 
https://petsc.org/main/manualpages/Mat/MatSetValuesLocal/

The local indices used for insertion here are indices relative to 
the block itself, and the L2G map for this matrix
has been rewritten to insert into that block in the larger matrix. 
Thus this looks like just calling MatSetValuesLocal()
on the smaller matrix block, but inserts correctly into the larger 
matrix.

Therefore, the code you write code in 2) could work equally well making the 
large matrix from 1), or independent smaller matrix blocks.

Does this make sense?

  Thanks,

 Matt


I will use MatSetLocalToGlobalMapping(K, row_mapping, column_mapping) to build 
the bottom left [P].





Q2) Can now I reset the mapping as MatSetLocalToGlobalMapping(K, 
column_mapping, row_mapping) to build the top right [P^T]?





Many thanks!



Kind regards,

Karthik.




















--
What most experimenters take for granted before they begin their experiments is 
infinitely more interesting than any results to which their experiments lead.
-- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/<http://www.cse.buffalo.edu/~knepley/>

Re: [petsc-users] question about MatSetLocalToGlobalMapping

[Karthikeyan Chockalingam - STFC UKRI via petsc-users]

Thank you, Matt. It took me a while but I understand how to work with 
submatrices.

(Below is simple stand-alone code which assembles two submatrices into a matrix 
of 5 x 5, which might be of use to someone else)

Q1) Looks like I don’t need to call MatrixAssembleBegin/End of on submatrices?
Q2) Though I have preallocated A for 5 x 5 – why do I still need to call

  MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE);
  else I get an error while assembling using MatSetValuesLocal()?
Q3) Am I calling MatRestoreLocalSubmatrix() at the right place and what does it 
do?
Q4) When assembling a large submatrix, can I call MatGetOwershipRange on the 
submatrix?

Best,
Karthik.





  #include 

  MatA;/* linear system matrix */

  PetscInt   i,j,m = 3,n = 3, P = 5, Q = 5;

  PetscErrorCode ierr;

  PetscScalarv;



  // create matrix A

  ierr = MatCreate(PETSC_COMM_WORLD,);CHKERRQ(ierr);

  ierr = MatSetType(A, MATMPIAIJ); CHKERRQ(ierr);

  ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,P,Q);CHKERRQ(ierr);

  ierr = MatSetFromOptions(A);CHKERRQ(ierr);

  ierr = MatMPIAIJSetPreallocation(A,5,NULL,5,NULL);CHKERRQ(ierr);

  //ierr = MatSeqAIJSetPreallocation(A,5,NULL);CHKERRQ(ierr);



  // local indices is always 0  1  2

 PetscInt  c = 3, global_col_indices[] = {0, 1, 2};

 PetscInt  r = 2, global_row_indices[] = {3, 4};



  ISLocalToGlobalMapping col_mapping, row_mapping;



  ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 1, c, global_col_indices, 
PETSC_COPY_VALUES, _mapping);

  ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 1, r, global_row_indices, 
PETSC_COPY_VALUES, _mapping);

  MatSetLocalToGlobalMapping(A,row_mapping, col_mapping);



  for (i=0; i
Date: Thursday, 20 April 2023 at 11:37
To: Chockalingam, Karthikeyan (STFC,DL,HC) 
Cc: petsc-users@mcs.anl.gov 
Subject: Re: [petsc-users] question about MatSetLocalToGlobalMapping
On Thu, Apr 20, 2023 at 6:13 AM Karthikeyan Chockalingam - STFC UKRI via 
petsc-users mailto:petsc-users@mcs.anl.gov>> wrote:
Hello,

I created a new thread, thought would it be more appropriate (and is a 
continuation of my previous post). I want to construct the below K matrix 
(which is composed of submatrices)

K = [A P^T
   P   0]

Where K is of type MatMPIAIJ. I first constructed the top left [A] using 
MatSetValues().

Now, I would like to construct the bottom left [p] and top right [p^T] using 
MatSetValuesLocal().

To use  MatSetValuesLocal(),  I first have to create a local-to-global mapping 
using ISLocalToGlobalMappingCreate. I have created two mapping row_mapping and 
column_mapping.

I do not understand why they are not the same map. Maybe I was unclear before. 
It looks like you have two fields, say phi and lambda, where lambda is a 
Lagrange multiplier imposing some constraint. Then you get a saddle point like 
this. You can imagine matrices

  (phi, phi)--> A
  (phi, lambda) --> P^T
  (lambda, phi) --> P

So you make a L2G map for the phi field and the lambda field. Oh, you are 
calling them row and col map, but they are my phi and lambda
maps. I do not like the row and col names since in P they reverse.

Q1) At what point should I declare MatSetLocalToGlobalMapping – is it just 
before I use MatSetValuesLocal()?

Okay, it is good you are asking this because my thinking was somewhat confused. 
I think the precise steps are:

  1) Create the large saddle point matrix K

1a) We must call 
https://petsc.org/main/manualpages/Mat/MatSetLocalToGlobalMapping/ on it. In 
the simplest case, this just maps
   the local rows numbers [0, Nrows) to the global rows numbers 
[rowStart, rowStart + Nrows).

  2) To form each piece:

2a) Extract that block using 
https://petsc.org/main/manualpages/Mat/MatGetLocalSubMatrix/

   This gives back a Mat object that you subsequently restore using 
https://petsc.org/main/manualpages/Mat/MatRestoreLocalSubMatrix/

 2b) Insert values using 
https://petsc.org/main/manualpages/Mat/MatSetValuesLocal/

The local indices used for insertion here are indices relative to 
the block itself, and the L2G map for this matrix
has been rewritten to insert into that block in the larger matrix. 
Thus this looks like just calling MatSetValuesLocal()
on the smaller matrix block, but inserts correctly into the larger 
matrix.

Therefore, the code you write code in 2) could work equally well making the 
large matrix from 1), or independent smaller matrix blocks.

Does this make sense?

  Thanks,

 Matt

I will use MatSetLocalToGlobalMapping(K, row_mapping, column_mapping) to build 
the bottom left [P].


Q2) Can now I reset the mapping as MatSetLocalToGlobalMapping(K, 
column_mapping, row_mapping) to build the top right [P^T]?


Many thanks!

Kind regards,
Karthik.











--
What most experimenters take for granted before they begin th

Re: [petsc-users] question about MatSetLocalToGlobalMapping

[Pierre Jolivet]

> On 20 Apr 2023, at 10:28 PM, Zhang, Hong  wrote:
> 
> Pierre,
> 1) Is there any hope to get PDIPDM to use a MatNest?
> 
> KKT matrix is indefinite and ill-conditioned, which must be solved using a 
> direct matrix factorization method.

But you are using PCBJACOBI in the paper you attached?
In any case, there are many such systems, e.g., a discretization of Stokes 
equations, that can be solved with something else than a direct factorization.

> For the current implementation, we use MUMPS Cholesky as default. To use 
> MatNest, what direct solver to use, SCHUR_FACTOR? I do not know how to get it 
> work. 

On the one hand, MatNest can efficiently convert to AIJ or SBAIJ if you want to 
stick to PCLU or PCCHOLESKY.
On the other hand, it allows to easily switch to PCFIELDSPLIT which can be used 
to solve saddle-point problems.

> 2) Is this fixed 
> https://lists.mcs.anl.gov/pipermail/petsc-dev/2020-September/026398.html ?
> I cannot get users to transition away from Ipopt because of these two missing 
> features.
> 
> The existing pdipm is the result of a MS student intern project. None of us 
> involved are experts on the optimization solvers. We made a straightforward 
> parallelization of Ipopt. It indeed needs further work, e.g., more features, 
> better matrix storage, convergence criteria... To our knowledge, parallel 
> pdipm is not available other than our pdipm.

Ipopt can use MUMPS and PARDISO internally, so it’s in some sense parallel 
(using shared memory).
Also, this is not a very potent selling point.
My users that are satisfied with Ipopt as a "non-parallel" black box don’t want 
to have to touch part of their code just to stick it in a parallel black box 
which is limited to the same kind of linear solver and which has severe 
limitations with respect to Hessian/Jacobian/constraint distributions.

Thanks,
Pierre

> We should improve our pdipm. 
> Hong
> 
>> On 20 Apr 2023, at 5:47 PM, Zhang, Hong via petsc-users 
>> mailto:petsc-users@mcs.anl.gov>> wrote:
>> 
>> Karthik,
>> We built a KKT matrix in  TaoSNESJacobian_PDIPM() (see 
>> petsc/src/tao/constrained/impls/ipm/pdipm.c) which assembles several small 
>> matrices into a large KKT matrix in mpiaij format. You could take the same 
>> approach to insert P and P^T into your K.
>> FYI, I attached our paper.
>> Hong
>>  
>> From: petsc-users > <mailto:petsc-users-boun...@mcs.anl.gov>> on behalf of Matthew Knepley 
>> mailto:knep...@gmail.com>>
>> Sent: Thursday, April 20, 2023 5:37 AM
>> To: Karthikeyan Chockalingam - STFC UKRI 
>> > <mailto:karthikeyan.chockalin...@stfc.ac.uk>>
>> Cc: petsc-users@mcs.anl.gov <mailto:petsc-users@mcs.anl.gov> 
>> mailto:petsc-users@mcs.anl.gov>>
>> Subject: Re: [petsc-users] question about MatSetLocalToGlobalMapping
>>  
>> On Thu, Apr 20, 2023 at 6:13 AM Karthikeyan Chockalingam - STFC UKRI via 
>> petsc-users mailto:petsc-users@mcs.anl.gov>> wrote:
>> Hello,
>>  
>> I created a new thread, thought would it be more appropriate (and is a 
>> continuation of my previous post). I want to construct the below K matrix 
>> (which is composed of submatrices)
>>  
>> K = [A P^T
>>P   0]
>>  
>> Where K is of type MatMPIAIJ. I first constructed the top left [A] using 
>> MatSetValues().
>>  
>> Now, I would like to construct the bottom left [p] and top right [p^T] using 
>> MatSetValuesLocal().
>>  
>> To use  MatSetValuesLocal(),  I first have to create a local-to-global 
>> mapping using ISLocalToGlobalMappingCreate. I have created two mapping 
>> row_mapping and column_mapping.
>> 
>> I do not understand why they are not the same map. Maybe I was unclear 
>> before. It looks like you have two fields, say phi and lambda, where lambda 
>> is a Lagrange multiplier imposing some constraint. Then you get a saddle 
>> point like this. You can imagine matrices
>> 
>>   (phi, phi)--> A
>>   (phi, lambda) --> P^T
>>   (lambda, phi) --> P
>> 
>> So you make a L2G map for the phi field and the lambda field. Oh, you are 
>> calling them row and col map, but they are my phi and lambda
>> maps. I do not like the row and col names since in P they reverse.
>>  
>> Q1) At what point should I declare MatSetLocalToGlobalMapping – is it just 
>> before I use MatSetValuesLocal()?
>> 
>> Okay, it is good you are asking this because my thinking was somewhat 
>> confused. I think the precise steps are:
>> 
>>   1) Create the large saddle point matrix K
>> 
>> 1a) We must call 
>> https://petsc.or

Re: [petsc-users] question about MatSetLocalToGlobalMapping

[Zhang, Hong via petsc-users]

Pierre,
1) Is there any hope to get PDIPDM to use a MatNest?

KKT matrix is indefinite and ill-conditioned, which must be solved using a 
direct matrix factorization method. For the current implementation, we use 
MUMPS Cholesky as default. To use MatNest, what direct solver to use, 
SCHUR_FACTOR? I do not know how to get it work.

2) Is this fixed 
https://lists.mcs.anl.gov/pipermail/petsc-dev/2020-September/026398.html ?
I cannot get users to transition away from Ipopt because of these two missing 
features.

The existing pdipm is the result of a MS student intern project. None of us 
involved are experts on the optimization solvers. We made a straightforward 
parallelization of Ipopt. It indeed needs further work, e.g., more features, 
better matrix storage, convergence criteria... To our knowledge, parallel pdipm 
is not available other than our pdipm. We should improve our pdipm.
Hong

On 20 Apr 2023, at 5:47 PM, Zhang, Hong via petsc-users 
 wrote:

Karthik,
We built a KKT matrix in  TaoSNESJacobian_PDIPM() (see 
petsc/src/tao/constrained/impls/ipm/pdipm.c) which assembles several small 
matrices into a large KKT matrix in mpiaij format. You could take the same 
approach to insert P and P^T into your K.
FYI, I attached our paper.
Hong

From: petsc-users  on behalf of Matthew 
Knepley 
Sent: Thursday, April 20, 2023 5:37 AM
To: Karthikeyan Chockalingam - STFC UKRI 
Cc: petsc-users@mcs.anl.gov 
Subject: Re: [petsc-users] question about MatSetLocalToGlobalMapping

On Thu, Apr 20, 2023 at 6:13 AM Karthikeyan Chockalingam - STFC UKRI via 
petsc-users mailto:petsc-users@mcs.anl.gov>> wrote:
Hello,

I created a new thread, thought would it be more appropriate (and is a 
continuation of my previous post). I want to construct the below K matrix 
(which is composed of submatrices)

K = [A P^T
   P   0]

Where K is of type MatMPIAIJ. I first constructed the top left [A] using 
MatSetValues().

Now, I would like to construct the bottom left [p] and top right [p^T] using 
MatSetValuesLocal().

To use  MatSetValuesLocal(),  I first have to create a local-to-global mapping 
using ISLocalToGlobalMappingCreate. I have created two mapping row_mapping and 
column_mapping.

I do not understand why they are not the same map. Maybe I was unclear before. 
It looks like you have two fields, say phi and lambda, where lambda is a 
Lagrange multiplier imposing some constraint. Then you get a saddle point like 
this. You can imagine matrices

  (phi, phi)--> A
  (phi, lambda) --> P^T
  (lambda, phi) --> P

So you make a L2G map for the phi field and the lambda field. Oh, you are 
calling them row and col map, but they are my phi and lambda
maps. I do not like the row and col names since in P they reverse.

Q1) At what point should I declare MatSetLocalToGlobalMapping – is it just 
before I use MatSetValuesLocal()?

Okay, it is good you are asking this because my thinking was somewhat confused. 
I think the precise steps are:

  1) Create the large saddle point matrix K

1a) We must call 
https://petsc.org/main/manualpages/Mat/MatSetLocalToGlobalMapping/ on it. In 
the simplest case, this just maps
   the local rows numbers [0, Nrows) to the global rows numbers 
[rowStart, rowStart + Nrows).

  2) To form each piece:

2a) Extract that block using 
https://petsc.org/main/manualpages/Mat/MatGetLocalSubMatrix/

   This gives back a Mat object that you subsequently restore using 
https://petsc.org/main/manualpages/Mat/MatRestoreLocalSubMatrix/

 2b) Insert values using 
https://petsc.org/main/manualpages/Mat/MatSetValuesLocal/

The local indices used for insertion here are indices relative to 
the block itself, and the L2G map for this matrix
has been rewritten to insert into that block in the larger matrix. 
Thus this looks like just calling MatSetValuesLocal()
on the smaller matrix block, but inserts correctly into the larger 
matrix.

Therefore, the code you write code in 2) could work equally well making the 
large matrix from 1), or independent smaller matrix blocks.

Does this make sense?

  Thanks,

 Matt

I will use MatSetLocalToGlobalMapping(K, row_mapping, column_mapping) to build 
the bottom left [P].


Q2) Can now I reset the mapping as MatSetLocalToGlobalMapping(K, 
column_mapping, row_mapping) to build the top right [P^T]?


Many thanks!

Kind regards,
Karthik.











--
What most experimenters take for granted before they begin their experiments is 
infinitely more interesting than any results to which their experiments lead.
-- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/<http://www.cse.buffalo.edu/~knepley/>

Re: [petsc-users] question about MatSetLocalToGlobalMapping

[Pierre Jolivet]

Hong,
1) Is there any hope to get PDIPDM to use a MatNest?
2) Is this fixed 
https://lists.mcs.anl.gov/pipermail/petsc-dev/2020-September/026398.html ?
I cannot get users to transition away from Ipopt because of these two missing 
features.

Thanks,
Pierre

> On 20 Apr 2023, at 5:47 PM, Zhang, Hong via petsc-users 
>  wrote:
> 
> Karthik,
> We built a KKT matrix in  TaoSNESJacobian_PDIPM() (see 
> petsc/src/tao/constrained/impls/ipm/pdipm.c) which assembles several small 
> matrices into a large KKT matrix in mpiaij format. You could take the same 
> approach to insert P and P^T into your K.
> FYI, I attached our paper.
> Hong
> From: petsc-users  on behalf of Matthew 
> Knepley 
> Sent: Thursday, April 20, 2023 5:37 AM
> To: Karthikeyan Chockalingam - STFC UKRI 
> Cc: petsc-users@mcs.anl.gov 
> Subject: Re: [petsc-users] question about MatSetLocalToGlobalMapping
>  
> On Thu, Apr 20, 2023 at 6:13 AM Karthikeyan Chockalingam - STFC UKRI via 
> petsc-users mailto:petsc-users@mcs.anl.gov>> wrote:
> Hello,
>  
> I created a new thread, thought would it be more appropriate (and is a 
> continuation of my previous post). I want to construct the below K matrix 
> (which is composed of submatrices)
>  
> K = [A P^T
>P   0]
>  
> Where K is of type MatMPIAIJ. I first constructed the top left [A] using 
> MatSetValues().
>  
> Now, I would like to construct the bottom left [p] and top right [p^T] using 
> MatSetValuesLocal().
>  
> To use  MatSetValuesLocal(),  I first have to create a local-to-global 
> mapping using ISLocalToGlobalMappingCreate. I have created two mapping 
> row_mapping and column_mapping.
> 
> I do not understand why they are not the same map. Maybe I was unclear 
> before. It looks like you have two fields, say phi and lambda, where lambda 
> is a Lagrange multiplier imposing some constraint. Then you get a saddle 
> point like this. You can imagine matrices
> 
>   (phi, phi)--> A
>   (phi, lambda) --> P^T
>   (lambda, phi) --> P
> 
> So you make a L2G map for the phi field and the lambda field. Oh, you are 
> calling them row and col map, but they are my phi and lambda
> maps. I do not like the row and col names since in P they reverse.
>  
> Q1) At what point should I declare MatSetLocalToGlobalMapping – is it just 
> before I use MatSetValuesLocal()?
> 
> Okay, it is good you are asking this because my thinking was somewhat 
> confused. I think the precise steps are:
> 
>   1) Create the large saddle point matrix K
> 
> 1a) We must call 
> https://petsc.org/main/manualpages/Mat/MatSetLocalToGlobalMapping/ on it. In 
> the simplest case, this just maps
>the local rows numbers [0, Nrows) to the global rows numbers 
> [rowStart, rowStart + Nrows).
> 
>   2) To form each piece:
> 
> 2a) Extract that block using 
> https://petsc.org/main/manualpages/Mat/MatGetLocalSubMatrix/
> 
>This gives back a Mat object that you subsequently restore using 
> https://petsc.org/main/manualpages/Mat/MatRestoreLocalSubMatrix/
> 
>  2b) Insert values using 
> https://petsc.org/main/manualpages/Mat/MatSetValuesLocal/
> 
> The local indices used for insertion here are indices relative to 
> the block itself, and the L2G map for this matrix
> has been rewritten to insert into that block in the larger 
> matrix. Thus this looks like just calling MatSetValuesLocal()
> on the smaller matrix block, but inserts correctly into the 
> larger matrix.
> 
> Therefore, the code you write code in 2) could work equally well making the 
> large matrix from 1), or independent smaller matrix blocks.
> 
> Does this make sense?
> 
>   Thanks,
> 
>  Matt
>  
> I will use MatSetLocalToGlobalMapping(K, row_mapping, column_mapping) to 
> build the bottom left [P].
>  
>  
> Q2) Can now I reset the mapping as MatSetLocalToGlobalMapping(K, 
> column_mapping, row_mapping) to build the top right [P^T]? 
>  
>  
> Many thanks!
>  
> Kind regards,
> Karthik.
>  
>  
>  
>  
>  
>  
>  
>  
>  
> 
> 
> -- 
> What most experimenters take for granted before they begin their experiments 
> is infinitely more interesting than any results to which their experiments 
> lead.
> -- Norbert Wiener
> 
> https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
>  interior point method for the solution of dynamic.pdf>

Re: [petsc-users] question about MatSetLocalToGlobalMapping

[Matthew Knepley]

On Thu, Apr 20, 2023 at 6:13 AM Karthikeyan Chockalingam - STFC UKRI via
petsc-users  wrote:

> Hello,
>
>
>
> I created a new thread, thought would it be more appropriate (and is a
> continuation of my previous post). I want to construct the below K matrix
> (which is composed of submatrices)
>
>
>
> K = [A P^T
>
>P   0]
>
>
>
> Where K is of type MatMPIAIJ. I first constructed the top left [A] using
> MatSetValues().
>
>
>
> Now, I would like to construct the bottom left [p] and top right [p^T] using
> MatSetValuesLocal().
>
>
>
> To use  MatSetValuesLocal(),  I first have to create a local-to-global
> mapping using ISLocalToGlobalMappingCreate. I have created two mapping
> row_mapping and column_mapping.
>

I do not understand why they are not the same map. Maybe I was unclear
before. It looks like you have two fields, say phi and lambda, where lambda
is a Lagrange multiplier imposing some constraint. Then you get a saddle
point like this. You can imagine matrices

  (phi, phi)--> A
  (phi, lambda) --> P^T
  (lambda, phi) --> P

So you make a L2G map for the phi field and the lambda field. Oh, you are
calling them row and col map, but they are my phi and lambda
maps. I do not like the row and col names since in P they reverse.


> Q1) At what point should I declare MatSetLocalToGlobalMapping – is it
> just before I use MatSetValuesLocal()?
>

Okay, it is good you are asking this because my thinking was somewhat
confused. I think the precise steps are:

  1) Create the large saddle point matrix K

1a) We must call
https://petsc.org/main/manualpages/Mat/MatSetLocalToGlobalMapping/ on it.
In the simplest case, this just maps
   the local rows numbers [0, Nrows) to the global rows numbers
[rowStart, rowStart + Nrows).

  2) To form each piece:

2a) Extract that block using
https://petsc.org/main/manualpages/Mat/MatGetLocalSubMatrix/

   This gives back a Mat object that you subsequently restore using
https://petsc.org/main/manualpages/Mat/MatRestoreLocalSubMatrix/

 2b) Insert values using
https://petsc.org/main/manualpages/Mat/MatSetValuesLocal/

The local indices used for insertion here are indices relative
to the block itself, and the L2G map for this matrix
has been rewritten to insert into that block in the larger
matrix. Thus this looks like just calling MatSetValuesLocal()
on the smaller matrix block, but inserts correctly into the
larger matrix.

Therefore, the code you write code in 2) could work equally well making the
large matrix from 1), or independent smaller matrix blocks.

Does this make sense?

  Thanks,

 Matt


> I will use MatSetLocalToGlobalMapping(K, row_mapping, column_mapping) to
> build the bottom left [P].
>
>
>
>
>
> Q2) Can now I reset the mapping as MatSetLocalToGlobalMapping(K,
> column_mapping, row_mapping) to build the top right [P^T]?
>
>
>
>
>
> Many thanks!
>
>
>
> Kind regards,
>
> Karthik.
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>


-- 
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/ 

[petsc-users] question about MatSetLocalToGlobalMapping

[Karthikeyan Chockalingam - STFC UKRI via petsc-users]

Hello,

I created a new thread, thought would it be more appropriate (and is a 
continuation of my previous post). I want to construct the below K matrix 
(which is composed of submatrices)

K = [A P^T
   P   0]

Where K is of type MatMPIAIJ. I first constructed the top left [A] using 
MatSetValues().

Now, I would like to construct the bottom left [p] and top right [p^T] using 
MatSetValuesLocal().

To use  MatSetValuesLocal(),  I first have to create a local-to-global mapping 
using ISLocalToGlobalMappingCreate. I have created two mapping row_mapping and 
column_mapping.

Q1) At what point should I declare MatSetLocalToGlobalMapping – is it just 
before I use MatSetValuesLocal()?


I will use MatSetLocalToGlobalMapping(K, row_mapping, column_mapping) to build 
the bottom left [P].


Q2) Can now I reset the mapping as MatSetLocalToGlobalMapping(K, 
column_mapping, row_mapping) to build the top right [P^T]?


Many thanks!

Kind regards,
Karthik.