Setting BC on a geometry inside another geometry

Hello community,

I have been trying to create a pipeline in which I load STL geometry and then create a box around that geometry so I can set boundary conditions to the sides such as in a virtual wind tunnel. I however need to set a no-slip boundary to both sides of the tunnel but as well to the geometry surface inside. I have created a geometry like this: )

blades = Tessellation.from_stl(geom_path + "/blades.stl", airtight=True)

    # normalize meshes
    def normalize_mesh(mesh, center, scale):
        mesh = mesh.translate([-c for c in center])
        # mesh = mesh.scale(scale) # not important right now
        return mesh

    idx = 0
    bnds = {}
    for key, value in blades.bounds.bound_ranges.items():
        bnds[blades.dims[idx]] = value
        idx += 1

    #  normalization
    def get_center(a, b):
        return ( a + b ) / 2

    x_offset = get_center(bnds['x'][0],bnds['x'][1])
    y_offset = get_center(bnds['y'][0],bnds['y'][1])
    z_offset = get_center(bnds['z'][0],bnds['z'][1])
    center = (x_offset, y_offset, z_offset)
    blades = normalize_mesh(blades, center, 1)
    # expecting that the geometry center is now set to (0, 0, 0)
    channel_length = (-400, 100)
    channel_width = (-250, 250)
    channel_height = (-250, 250)
    box_bounds = {x: channel_length, y: channel_width, z: channel_height}

    # define interior geometry, without blades
    rec = Box(
        (channel_width[0], channel_length[0], channel_height[0]),
        (channel_width[1], channel_length[1], channel_height[1]),
        parameterization=OrderedParameterization(time_range, key=t_symbol)
    ) - blades

    geo = rec + blades

On this, I have created the no slip for walls like this and it works:

    noslipBC = PointwiseBoundaryConstraint(
        nodes=nodes,
        geometry=geo,
        outvar={"u": 0, "v": 0, "w": 0, "p" : 0},
        batch_size=cfg.batch_size.initial_condition,
        #bounds=box_bounds,
        lambda_weighting={"u": 100, "v": 100, "w": 100, "p": 100},
        # criteria for all side walls
        criteria=And((y > channel_length[0]), 
                    (y < channel_length[1]), 
                    Or(
                        Or(Eq(x, channel_width[0]), Eq(x, channel_width[1])), 
                        Or(Eq(z, channel_height[0]), Eq(z, channel_height[1]))
                    )
                ),
    )
    ic_domain.add_constraint(noslipBC, "noslipBC")
    window_domain.add_constraint(noslipBC, "noslipBC")

But if I want to select the partial geometry, such as here:

    bladesBC = PointwiseBoundaryConstraint(
        nodes=nodes,
        geometry=blades,
        outvar={"u": 0, "v": 0, "w": 0},
        batch_size=cfg.batch_size.initial_condition,
        lambda_weighting={"u": 100, "v": 100, "w": 100, "p": 100},
    )
    ic_domain.add_constraint(bladesBC, "bladesBC")
    window_domain.add_constraint(bladesBC, "bladesBC")

The whole scripts ends up with failed unrolling graph:

####################################
could not unroll graph!
This is probably because you are asking to compute a value that is not an output of any node
####################################
invar: [x, y, z, sdf, area]
requested var: [u, v, w]
computable var: [x, y, z, sdf, area, continuity]
####################################
Nodes in graph: 
node: Sympy Node: continuity
evaluate: SympyToTorch
inputs: []
derivatives: [u__x, v__y, w__z]
outputs: [continuity]
optimize: False
node: Sympy Node: momentum_x
evaluate: SympyToTorch
inputs: [u, v, w]
derivatives: [p__x, u__t, u__x, u__x__x, u__y, u__y__y, u__z, u__z__z]
outputs: [momentum_x]
optimize: False
node: Sympy Node: momentum_y
evaluate: SympyToTorch
inputs: [u, v, w]
derivatives: [p__y, v__t, v__x, v__x__x, v__y, v__y__y, v__z, v__z__z]
outputs: [momentum_y]
optimize: False
node: Sympy Node: momentum_z
evaluate: SympyToTorch
inputs: [u, v, w]
derivatives: [p__z, w__t, w__x, w__x__x, w__y, w__y__y, w__z, w__z__z]
outputs: [momentum_z]
optimize: False
node: Arch Node: time_window_network
evaluate: MovingTimeWindowArch
inputs: [x, y, z, t]
derivatives: []
outputs: [u, v, w, p, u_prev_step, v_prev_step, w_prev_step, p_prev_step, u_prev_step_diff, v_prev_step_diff, w_prev_step_diff, p_prev_step_diff]
optimize: True
####################################
Error executing job with overrides: []
Traceback (most recent call last):
  File "/root/Desktop/workspace/examples/taylor_green/turbine.py", line 194, in run
    bladesBC = PointwiseInteriorConstraint(
  File "/opt/conda/lib/python3.8/site-packages/modulus-22.9-py3.8.egg/modulus/domain/constraint/continuous.py", line 493, in __init__
    super().__init__(
  File "/opt/conda/lib/python3.8/site-packages/modulus-22.9-py3.8.egg/modulus/domain/constraint/constraint.py", line 55, in __init__
    self.model = Graph(
  File "/opt/conda/lib/python3.8/site-packages/modulus-22.9-py3.8.egg/modulus/graph.py", line 94, in __init__
    raise RuntimeError("Failed Unrolling Graph")
RuntimeError: Failed Unrolling Graph

Don’t mind the root dir, it is docker container :)

The question here is, how am I supposed to add boundary condition for the geometry I have placed inside the box? I have also similar problem if I use rec geometry alone, the I get - loss went to nans. Adding code here as well:

    interior = PointwiseInteriorConstraint(
        nodes=nodes,
        geometry=rec,
        outvar={"continuity": 0, "momentum_x": 0, "momentum_y": 0, "momentum_z": 0},
        bounds=box_bounds,
        batch_size=4000,
        loss=CausalLossNorm(eps=1.0),
        fixed_dataset=False,
        shuffle=False
    )
    ic_domain.add_constraint(interior, name="interior")
    window_domain.add_constraint(interior, name="interior")

I guess I just don’t fully understand the logic behing geometry in modulus. Thanks for any insigths :)

Well so I have found a reason after all and here is what has caused the weird behaviour:

I haven’t noticed, that if I use parametrization for one of my geometries, I should have do it for all of the geometries. So simply replacing this line:

blades = Tessellation.from_stl(geom_path + "/blades.stl", airtight=True)

by this:

blades = Tessellation.from_stl(geom_path + "/blades.stl", airtight=True,
        parameterization=OrderedParameterization(time_range, key=t_symbol))

has solved the whole problem.

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.