Thanks, Matt. I am using the following steps to build a local to global
mapping.
Step 1) PetscSectionCreate ();
PetscSectionSetNumFields();
PetscSectionSetChart ();
//Set dof for each node
PetscSectionSetup (s);
Step 2)
PetscCall(DMGetLocalToGlobalMapping(dm, <ogm));
PetscCall(ISLocalToGlobalMappingGetIndices(ltogm, &g_idx));
For 2D mesh from the example,
https://wg-beginners.readthedocs.io/en/latest/tutorials/dm/plex/introductory_tutorial_plex.html
the map from rank 0 (MPI) and rank 1(MPI) matches well with the global node
ordering.
But When I tried a 3D gmsh, (2 hexahedra (12 nodes ) were split into 12
tetrahedra). ( 2 processors in total ) Each processor handles 9 nodes
separately. When I checked the local to global mapping, It seems the
mapping is not right.
For example, the mapping array is
Local Global (left: local node number; right: global node number).
provided.
0 0
1 3
2 5
3 1
4 6
5 8
6 2
7 9
8 11
But the coordinates between the local (left) and global nodes (right) are
not the same, meaning they are not matching. Did I miss something? In
addition, the simple 3D gmsh file has been attached if you need.
*Gmsh file:*
$MeshFormat
4.1 0 8
$EndMeshFormat
$PhysicalNames
3
2 27 "Excit"
2 28 "PEC"
3 29 "volume"
$EndPhysicalNames
$Entities
8 12 6 1
1 0 0 0 0
2 0 0.25 0 0
3 0 0.25 0.25 0
4 0 0 0.25 0
5 1 0 0 0
6 1 0.25 0 0
10 1 0.25 0.25 0
14 1 0 0.25 0
1 0 0 0 0 0.25 0 0 2 1 -2
2 0 0.25 0 0 0.25 0.25 0 2 2 -3
3 0 0 0.25 0 0.25 0.25 0 2 3 -4
4 0 0 0 0 0 0.25 0 2 4 -1
6 1 0 0 1 0.25 0 0 2 5 -6
7 1 0.25 0 1 0.25 0.25 0 2 6 -10
8 1 0 0.25 1 0.25 0.25 0 2 10 -14
9 1 0 0 1 0 0.25 0 2 14 -5
11 0 0 0 1 0 0 0 2 1 -5
12 0 0.25 0 1 0.25 0 0 2 2 -6
16 0 0.25 0.25 1 0.25 0.25 0 2 3 -10
20 0 0 0.25 1 0 0.25 0 2 4 -14
1 0 0 0 0 0.25 0.25 1 27 4 1 2 3 4
13 0 0 0 1 0.25 0 1 28 4 1 12 -6 -11
17 0 0.25 0 1 0.25 0.25 0 4 2 16 -7 -12
21 0 0 0.25 1 0.25 0.25 1 28 4 3 20 -8 -16
25 0 0 0 1 0 0.25 0 4 4 11 -9 -20
26 1 0 0 1 0.25 0.25 1 28 4 6 7 8 9
1 0 0 0 1 0.25 0.25 1 29 6 -1 26 13 17 21 25
$EndEntities
$Nodes
17 12 1 12
0 1 0 1
1
0 0 0
0 2 0 1
2
0 0.25 0
0 3 0 1
3
0 0.25 0.25
0 4 0 1
4
0 0 0.25
0 5 0 1
5
1 0 0
0 6 0 1
6
1 0.25 0
0 10 0 1
7
1 0.25 0.25
0 14 0 1
8
1 0 0.25
1 11 0 1
9
0.5 0 0
1 12 0 1
10
0.5 0.25 0
1 16 0 1
11
0.5 0.25 0.25
1 20 0 1
12
0.5 0 0.25
2 1 0 0
2 13 0 0
2 21 0 0
2 26 0 0
3 1 0 0
$EndNodes
$Elements
5 24 1 24
2 1 2 2
1 1 2 4
2 4 2 3
2 13 2 4
3 9 1 2
4 9 2 10
5 5 9 10
6 5 10 6
2 21 2 4
7 11 3 4
8 11 4 12
9 7 11 12
10 7 12 8
2 26 2 2
11 5 6 8
12 8 6 7
3 1 4 12
13 1 2 4 12
14 10 9 12 2
15 2 9 12 1
16 9 10 12 8
17 6 5 8 10
18 10 5 8 9
19 4 2 3 11
20 10 12 11 2
21 2 12 11 4
22 12 10 11 7
23 6 8 7 10
24 10 8 7 12
$EndElements
On Wed, May 17, 2023 at 7:30 PM Matthew Knepley <knep...@gmail.com> wrote:
> On Wed, May 17, 2023 at 6:58 PM neil liu <liufi...@gmail.com> wrote:
>
>> Dear Petsc developers,
>>
>> I am writing my own code to calculate the FEM matrix. The following is my
>> general framework,
>>
>> DMPlexCreateGmsh();
>> MPI_Comm_rank (Petsc_comm_world, &rank);
>> DMPlexDistribute (.., .., &dmDist);
>>
>> dm = dmDist;
>> //This can create separate dm s for different processors. (reordering.)
>>
>> MatCreate (Petsc_comm_world, &A)
>> // Loop over every tetrahedral element to calculate the local matrix for
>> each processor. Then we can get a local matrix A for each processor.
>>
>> *My question is : it seems we should build a global matrix B (assemble
>> all the As for each partition) and then transfer B to KSP. KSP will do the
>> parallelization correctly, right? *
>>
>
> I would not suggest this. The more common strategy is to assemble each
> element matrix directly into the
> global matrix B, by mapping the cell indices directly to global indices
> (rather than to local indices in the matrix A). You can do this in two
> stages. You can create a LocalToGlobalMapping in PETSc that maps
> every local index to a global index. Then you can assemble into B exactly
> as you would assemble into A by calling MatSetValuesLocal().
>
> DMPlex handles these mappings for you automatically, but I realize that it
> is a large number of things to buy into.
>
> Thanks,
>
> Matt
>
>
>> If that is right, I should define a whole domain matrix B before the
>> partitioning (MatCreate (Petsc_comm_world, &B); ), and then use
>> localtoglobal (which petsc function should I use? Do you have any
>> examples.) map to add A to B at the right positions (MatSetValues) ?
>>
>> Does that make sense?
>>
>> Thanks,
>>
>> Xiaodong
>>
>>
>>
>>
>>
>>
>>
>>
>
> --
> 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/>
>
