Interpolation between meshes

One of the key features of TMR is that it can be used to generate a hierarchy of meshes for multigrid. TMR uses an interpolation code from TACS to define interpolations between different mesh levels. This example demonstrates this capability and the ability to create higher-order meshes. Note that the meshes need not be uniformly refined/coarsened.

This code consists of the following steps:

• The geometry is created by importing a STEP file

• Create the fine octree mesh with 4 nodes per edge with the Gauss-Lobatto nodes

• Set the refinement array, refining the first element 4 times

• Balance the mesh and repartition and create the nodes on the fine mesh

• Make the coarse mesh and perform a random refinement

• Create variable map objects in TACS needed for the interpolation object

• Create the interpolations between the coarse and fine meshes

import os
import numpy as np
from mpi4py import MPI
from tmr import TMR
from tacs import TACS

# Set the communicator
comm = MPI.COMM_WORLD

# The fine octree forest
fine = None

stepfile = 'beam.stp'
if os.path.isfile(stepfile):
    # Load the geometry model
    geo = TMR.LoadModel(stepfile)

    # Mark the boundary condition faces
    faces = geo.getFaces()
    volumes = geo.getVolumes()
    faces[4].setSource(volumes[0], faces[5])

    # Create the mesh
    mesh = TMR.Mesh(comm, geo)

    # Set the meshing options
    opts = TMR.MeshOptions()
    opts.frontal_quality_factor = 1.25
    opts.num_smoothing_steps = 10
    opts.write_mesh_quality_histogram = 0

    # Create the surface mesh
    htarget = 4.0
    mesh.mesh(htarget, opts)

    # Create a model from the mesh
    model = mesh.createModelFromMesh()

    # Create the corresponding mesh topology from the mesh-model
    topo = TMR.Topology(comm, model)
    fine = TMR.OctForest(comm)
    fine.setTopology(topo)
else:
    conn = np.array([[0, 1, 3, 4, 6, 7, 9, 10],
                    [8, 11, 2, 5, 7, 10, 1, 4]], dtype=np.intc)
    fine = TMR.OctForest(comm)
    fine.setConnectivity(conn)

# Create the fine mesh
fine.createTrees(0)
refine = np.zeros(len(fine.getOctants()), dtype=np.intc)
if comm.rank == 0:
    refine[0] = 4
fine.refine(refine)
fine.balance(1)
fine.repartition()
fine.setMeshOrder(4, TMR.GAUSS_LOBATTO_POINTS)
fine.createNodes()

# Create a refinement array
octants = fine.getOctants()
refine = np.zeros(len(octants), dtype=np.intc)
for i in range(len(octants)):
    if i % 7 == 0:
        refine[i] = 2
    elif i % 5 == 0:
        refine[i] = -1

# Make the coarse tree finer than the fine tree for testing purposes
coarse = fine.duplicate()
coarse.refine(refine)
coarse.balance(0)
coarse.repartition()
coarse.setMeshOrder(2, TMR.GAUSS_LOBATTO_POINTS)
coarse.createNodes()

coarse_range = coarse.getNodeRange()
nc = coarse_range[comm.rank+1] - coarse_range[comm.rank]
coarse_map = TACS.VarMap(comm, nc)

fine_range = fine.getNodeRange()
nf = fine_range[comm.rank+1] - fine_range[comm.rank]
fine_map = TACS.VarMap(comm, nf)

# Create the two interpolations fine -> coarse and coarse -> fine
interp = TACS.VecInterp(coarse_map, fine_map, 1)
fine.createInterpolation(coarse, interp)
interp.initialize()

interp2 = TACS.VecInterp(fine_map, coarse_map, 1)
coarse.createInterpolation(fine, interp2)
interp2.initialize()