A Demo

Here, we show a demo for transferring solutions between two meshes of the unit square. We let the blue mesh participant run in parallel with two cores, while the green side is treated as a serial mesh.

import numpy as np
from mpi4py import MPI
from parpydtk2 import *

comm = MPI.COMM_WORLD

blue, green = create_imeshdb_pair(comm)

assert comm.size == 2
rank = comm.rank

For demo purpose, we construct the meshes globally.

# create blue meshes on all processes
cob = np.empty(shape=(16, 3), dtype=float, order='C')
dh = 0.33333333333333333
index = 0
x = 0.0
for i in range(4):
    y = 0.0
    for j in range(4):
        cob[index] = np.asarray([x, y, 0.0])
        index += 1
        y += dh
    x += dh

# global IDs, one based
bgids = np.arange(16, dtype='int32')
bgids += 1

The blue side has 16 nodes, the following is for the green side:

# create green meshes on all processes
cog = np.empty(shape=(36, 3), dtype=float, order='C')
dh = 0.2
index = 0
x = 0.0
for i in range(6):
    y = 0.0
    for j in range(6):
        cog[index] = np.asarray([x, y, 0.0])
        index += 1
        y += dh
    x += dh

# global IDs
ggids = np.arange(36, dtype='int32')
ggids += 1

The green participant has 36 nodes. The next step is to put the data in the two mesh databases.

# creating the mesh database
blue.begin_create()
# local vertices and global IDs
lcob = cob[8 * rank:8 * (rank + 1), :].copy()
lbgids = bgids[8 * rank:8 * (rank + 1)].copy()
blue.create_vertices(lcob)
blue.assign_gids(lbgids)
# do not use trivial global ID strategy
blue.finish_create(False)
blue.create_field('b')

As we can see, we equally distributed the mesh into the two cores as well as the corresponding global IDs. In line 54, the False flag indicates that the mesh database should use the user-provided global IDs.

Warning

Creating vertices and assigning global IDs must be called between begin_create() and finish_create()! Otherwise, exceptions are thrown.

Warning

Creating fields must be done after finish_create()!

Here is the treatment for the “serial” participant:

# NOTE that green is assumed to be serial mesh
green.begin_create()
# only create on master rank
if not rank:
    green.create_vertices(cog)
# since green is serial, we just use the trivial global IDs
green.finish_create()  # empty partition is resolve here
green.create_field('g')  # must after finish create

assert green.has_empty()

As we can see, only the master process has data.

Note

The duplicated node is handled inside finish_create()

With the two participants ready, we can now create our Mapper.

# create our analytical model, i.e. 10+sin(x)*cos(y)
bf = 10.0 + np.sin(lcob[:, 0]) * np.cos(lcob[:, 1])
gf = np.sin(cog[:, 0]) * np.cos(cog[:, 1]) + 10.0

# Construct our mapper
mapper = Mapper(blue=blue, green=green)

# using mmls, blue radius 1.0 green radius 0.6
mapper.method = MMLS
mapper.radius_b = 1.0
mapper.radius_g = 0.6
mapper.dimension = 2

# Mapper initialization region
mapper.begin_initialization()
mapper.register_coupling_fields(bf='b', gf='g', direct=B2G)
mapper.register_coupling_fields(bf='b', gf='g', direct=G2B)
mapper.end_initialization()

Line 69-70 just create an analytic model for error analysis. Line 77-80 are for parameters, for this case, we use MMLS (default) with blue side radius 1.0 and green side radius 0.6 for searching.

The important part is from line 83 to 86. Particularly speaking, the function register_coupling_fields(). It takes three parameters, where the first two are string tokens that represents the data fields in blue and green. The direct is to indicate the transfer direction, e.g. B2G stands for blue to green.

Warning

register_coupling_fields() must be called within begin_initialization() and end_initialization().

# NOTE that the following only runs on green master mesh
if not rank:
    green.assign_field('g', gf)
# Since green is serial and has empty partition, we must call this to
# resolve asynchronous values
green.resolve_empty_partitions('g')

The above section is to assign values on the green participant. Notice that it is a “serial” mesh, so we only assign values on the master process. But resolve the duplicated node is needed, this is done in line 93.

Finally, the solution transfer part is pretty straightforward:

# solution transfer region
mapper.begin_transfer()
mapper.transfer_data(bf='b', gf='g', direct=G2B)
err_b = (bf - blue.extract_field('b'))/10
mapper.transfer_data(bf='b', gf='g', direct=B2G)
err_g = (gf - green.extract_field('g'))/10
mapper.end_transfer()

comm.barrier()

print(rank, 'blue L2-error=%.3e' % (np.linalg.norm(err_b)/np.sqrt(err_b.size)))
if rank == 0:
    print(0, 'green L2-error=%.3e' % (np.linalg.norm(err_g)/np.sqrt(err_g.size)))

This code can be obtained here parallel2serial.py.