Skip to content

Commit b76d8fd

Browse files
fluidnumericsJoeclaude
fluidnumericsJoe
and
claude
authored
Variable sound speed in LinearEuler2D + HOHQMesh ISM/ISM-MM reader (#123)
* Make sound speed a per-node solution variable in LinearEuler2D Move c from a scalar model attribute to solution variable index 5 so the sound speed can vary in space. Volume flux, Riemann flux, and source for the new variable are identically zero, so c is held fixed in time at each location. - _t module: nvar=5, c metadata, drop %c scalar, flux/riemann/entropy/ no-normal-flow BC updated; SphericalSoundWave now takes c as an arg. - GPU Fortran: drop c parameter from interfaces and call sites; radiation BC wrapper now also passes the interior boundary buffer. - GPU C++ kernels: read c from solution[idof+4*ndof] (or boundary state), use max(cL,cR) for Lax-Friedrichs dissipation, emit zero flux for var 5; no-normal-flow copies c through; radiation preserves c from interior. - Examples: planewave_{propagation,reflection} carry their own c attribute (parent scalar is gone) and write solution(...,5) in both initial conditions and prescribed BCs; spherical_soundwave passes c into SphericalSoundWave. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * Add HOHQMesh ISM and ISM-MM reader for Mesh2D, with material tracking Adds `Read_HOHQMesh` to `Mesh2D_t` that auto-detects ISM vs ISM-MM from the first line of the file, parses nodes and per-element blocks, and reconstructs SELF's full mesh state — including side connectivity, which is not present in either format and must be derived from corner-node-id matching. - New fields on `Mesh2D_t`: `nMaterials`, `elemMaterial(:)`, `materialNames(:)`. Default-allocated by `Init` to a single "default" material so HOPr/structured readers stay unchanged. - HOHQMesh curve samples are at Chebyshev-Gauss-Lobatto nodes; the reader sets `quadrature = CHEBYSHEV_GAUSS_LOBATTO` and uses SELF's `ChebyshevQuadrature` for the TFI parametric coords. - Curved-edge interiors (polyOrder > 1) are filled via transfinite (Coons-patch) interpolation. Sample order is reversed on the North/West sides to match SELF's (ξ, η) convention. - Side connectivity rebuilt by hashing sorted corner-id pairs: fills `sideInfo(3:4,...)`, deduplicates global side ids, and marks boundary edges via a bc-name lookup table; `RecalculateFlip` resolves orientation. - Adds tests `mesh2d_read_ism` and `mesh2d_read_ismmm` plus the BoneAndMarrow ISM-MM fixture (and its .control file). The Block2D fixture is reused. Side accounting and material tables are validated against expected values. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * Add multi-material BoneAndMarrow LinearEuler2D example, fix TFI orientation New CPU example demonstrates the variable-sound-speed LinearEuler2D model on the ISM-MM BoneAndMarrow mesh. A small Gaussian pressure and density bump is seeded in the Muscle annulus, c is initialised per element from the parsed material name (Muscle=1.0, Bone=2.3, Marrow=0.92 in normalized units), and a CPU radiation BC handler is registered on top of the base class. All four outer-boundary sides remap to SELF_BC_RADIATION via ResetBoundaryConditionType. While bringing the example up I found the transfinite interpolation was reversing curve samples for the North and West sides. HOHQMesh's edgeMap (1,2;2,3;4,3;1,4) already runs sides 3 and 4 in SELF's +ξ and +η directions, so my reversal was folding curved elements at material interfaces (62 of 973 elements had Jacobians that changed sign). With the reversal removed, every Jacobian is positive and the simulation runs stably: entropy decays monotonically from 2.83e-7 to 1.97e-7 over 500 RK3 steps, c is preserved bitwise at the material values, and the acoustic pulse spreads cleanly through the heterogeneous medium and out through the radiation boundary. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * Mirror material-table alloc/free in GPU Mesh2D Init/Free The GPU Mesh2D duplicates host-array allocations inline rather than delegating to the parent Init/Free, so the new elemMaterial and materialNames fields would be left unallocated on a GPU build — Read_HOHQMesh's deallocate / assignment would then segfault. Replicate the parent's allocation pattern: allocate a single "default" material entry covering every element in Init, and deallocate the pair in Free. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * Apply fprettify formatting; unpin x86 CPU CI from franklin Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * Replace LLF with impedance-matched Riemann flux for LinearEuler2D The 2-D linear-Euler model with variable sound speed (variable 5) was using LLF / Rusanov with cmax = max(c_L, c_R). On a curved multi-material mesh with discontinuous c across element interfaces this is non-conservative for the pressure equation (each side's flux uses its own c^2; the central average mismatches both) and at high polynomial order the resulting aliasing instability blows the solution up the first time an acoustic wave encounters the c jump. Replaced with the exact (characteristic-decomposition) Riemann flux for linear acoustics with variable c. The normal-flux Jacobian eigenstructure is +c : right-going acoustic W_+ = rho0*c*u_n + p -c : left-going acoustic W_- = -rho0*c*u_n + p 0 : entropy density W_0 = rho' - p/c^2 0 : tangential vorticity u_t Upwinding each mode by its characteristic direction across the face (W_+|* = W_+|_L, W_-|* = W_-|_R) yields the impedance-matched interface state u_n* = (Z_L u_n,L + Z_R u_n,R + (p_L - p_R)) / (Z_L + Z_R) p* = (Z_R p_L + Z_L p_R + Z_L Z_R (u_n,L - u_n,R)) / (Z_L + Z_R) with Z = rho0*c. This is exact upwind / Godunov for the linearised acoustic system and reduces correctly to Fresnel-style reflection / transmission at a c-discontinuity. The pressure flux uses an averaged c^2 as a pragmatic treatment of the non-conservative product rho0*c^2*div(v) at the face; a fully path-conservative (Castro-Pares) treatment would use each side's own c^2 in its surface integral. Verified on the BoneAndMarrow example with N=7 controlDegree on GAUSS quadrature, dt=1e-4, endtime=40. The previous LLF scheme diverged at t~4 (the first wave-bone encounter) at N=5 and N=7 regardless of dt; the new flux runs cleanly through 40 sim seconds with pressure remaining bounded and entropy monotonically decreasing as the disturbance radiates out of the domain. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> * Inherit AdditionalInit in BoneAndMarrow example so GPU radiation BC actually runs The example overrode AdditionalInit and chained to AdditionalInit_LinearEuler2D_t (the CPU base). On a GPU build that bypassed AdditionalInit_LinearEuler2D (the GPU parent), so the GPU radiation BC wrapper hbc2d_Radiation_LinearEuler2D_GPU_wrapper was never registered. The example's own CPU radiation handler ran instead, modifying host buffers that the device pipeline never reads — leaving solution%extBoundary_gpu at outer boundary faces with c = 0. With c_R = 0 the impedance-matched Riemann flux degenerates: u_n* picks up the entire interior pressure (p_L / Z_L), p* collapses to zero. That's not a radiation condition — it's effectively a soft-wall reflector that inverts the wave amplitude. Visible in slices as strong reflections off the outer disk boundary. Fix is to remove the AdditionalInit override (and the now-redundant CPU radiation handler hbc2d_Radiation_lineareuler2d_boneandmarrow). The example type inherits AdditionalInit from lineareuler2d, which is the GPU parent's AdditionalInit_LinearEuler2D. That registers both the GPU no_normal_flow wrapper and the GPU radiation wrapper, and the wrapper correctly preserves c on the device-side extBoundary buffer. On the 4-radius disk with the Gaussian pulse, the bounded amplitude after 40 simulation seconds drops by a factor of ~3000 in peak |P| vs the broken-BC run, and visible boundary reflections vanish. The viz mp4 size shrinks from 1.3 MB to 477 KB — most of the wave content is genuinely radiating out instead of bouncing back. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com> --------- Co-authored-by: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
1 parent f5110f5 commit b76d8fd

16 files changed

Lines changed: 6226 additions & 89 deletions

‎.buildkite/x86-cpu-debug-tests.yml‎

Lines changed: 0 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -39,6 +39,5 @@ steps:
3939
slurm_ntasks: 2
4040
slurm_cpus_per_task: 8
4141
slurm_container_image: "docker://higherordermethods/selfish:latest-x86-none"
42-
slurm_nodelist: "franklin"
4342
agents:
4443
queue: "galapagos"

‎examples/CMakeLists.txt‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -65,6 +65,7 @@ add_fortran_tests (
6565
"linear_euler2d_spherical_soundwave_closeddomain.f90"
6666
"linear_euler2d_planewave_reflection.f90"
6767
"linear_euler2d_planewave_propagation.f90"
68+
"linear_euler2d_boneandmarrow.f90"
6869
"linear_shallow_water2d_nonormalflow.f90"
6970
"linear_shallow_water2d_planetaryrossby_wave.f90"
7071
"linear_shallow_water2d_kelvinwaves.f90"

‎examples/linear_euler2d_boneandmarrow.f90‎

Lines changed: 184 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,184 @@
1+
! //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// !
2+
!
3+
! Maintainers : support@fluidnumerics.com
4+
! Official Repository : https://github.com/FluidNumerics/self/
5+
!
6+
! Copyright © 2026 Fluid Numerics LLC
7+
!
8+
! Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
9+
!
10+
! 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
11+
!
12+
! 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in
13+
! the documentation and/or other materials provided with the distribution.
14+
!
15+
! 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from
16+
! this software without specific prior written permission.
17+
!
18+
! //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// !
19+
module lineareuler2d_boneandmarrow_model
20+
!! Heterogeneous-medium linear-acoustics test on the
21+
!! share/mesh/MultiMaterial2D/BoneAndMarrow.mesh ISM-MM mesh.
22+
!!
23+
!! The mesh tags every element with one of three materials:
24+
!! "Muscle" (background annulus, r in roughly [6.5, 14]),
25+
!! "Bone" (disk of radius 6.5 about (-0.5, 0)), and
26+
!! "Marrow" (smaller disk of radius 1.5 about (-1.5, 0.5))
27+
!! nested inside the bone region. We map each material to a
28+
!! representative sound speed and write that into solution(...,5),
29+
!! which is held fixed in time by the LinearEuler2D model.
30+
!!
31+
!! Initial condition: a small Gaussian pressure / density bump is
32+
!! placed in the Muscle region, well outside the bone region. The
33+
!! transient acoustic pulse propagates outward, refracts/reflects
34+
!! at the material interfaces (because c is discontinuous across
35+
!! them), and radiates out through the outer boundary.
36+
37+
use self_lineareuler2d
38+
39+
implicit none
40+
41+
type,extends(lineareuler2d) :: lineareuler2d_boneandmarrow
42+
! Sound speeds in normalized units, preserving real-tissue ratios
43+
! (bone is ~2.3x faster than muscle, marrow is ~0.92x).
44+
real(prec) :: c_muscle = 1.0_prec
45+
real(prec) :: c_bone = 2.3_prec
46+
real(prec) :: c_marrow = 0.92_prec
47+
real(prec) :: bump_x0 = 10.0_prec ! Pulse center, well into the muscle annulus
48+
real(prec) :: bump_y0 = 0.0_prec
49+
real(prec) :: bump_L = 0.6_prec ! Halfwidth (e-folding length) of the bump
50+
real(prec) :: bump_amp = 1.0e-3_prec ! Density-perturbation amplitude
51+
52+
contains
53+
54+
! AdditionalInit is inherited from `lineareuler2d`, which registers both
55+
! the no-normal-flow and radiation BCs using the correct CPU/GPU wrappers
56+
! for the current build. The radiation wrapper zeros (rho', u, v, p) and
57+
! preserves c on the same host/device buffer the rest of the pipeline
58+
! reads from, which is what this disk-shaped scatterer test requires.
59+
procedure :: setInitialCondition
60+
61+
endtype lineareuler2d_boneandmarrow
62+
63+
contains
64+
65+
subroutine setInitialCondition(this)
66+
!! Material-aware initial condition: stamp `c` from the per-
67+
!! element material id, then add a Gaussian rho/p bump only in
68+
!! Muscle elements (so the pulse starts cleanly in a uniform
69+
!! medium and the bone/marrow inclusions are seen as scatterers).
70+
implicit none
71+
class(lineareuler2d_boneandmarrow),intent(inout) :: this
72+
integer :: i,j,iel,matid
73+
real(prec) :: c_mat,x,y,r2,shape
74+
character(LEN=64) :: matname
75+
76+
do iel = 1,this%mesh%nElem
77+
matid = this%mesh%elemMaterial(iel)
78+
matname = this%mesh%materialNames(matid)
79+
select case(trim(matname))
80+
case("Muscle"); c_mat = this%c_muscle
81+
case("Bone"); c_mat = this%c_bone
82+
case("Marrow"); c_mat = this%c_marrow
83+
case default; c_mat = this%c_muscle
84+
endselect
85+
86+
do j = 1,this%solution%N+1
87+
do i = 1,this%solution%N+1
88+
x = this%geometry%x%interior(i,j,iel,1,1)
89+
y = this%geometry%x%interior(i,j,iel,1,2)
90+
91+
if(trim(matname) == "Muscle") then
92+
r2 = (x-this%bump_x0)**2+(y-this%bump_y0)**2
93+
shape = this%bump_amp*exp(-r2/(this%bump_L*this%bump_L))
94+
else
95+
shape = 0.0_prec
96+
endif
97+
98+
this%solution%interior(i,j,iel,1) = shape ! density perturbation
99+
this%solution%interior(i,j,iel,2) = 0.0_prec ! u
100+
this%solution%interior(i,j,iel,3) = 0.0_prec ! v
101+
this%solution%interior(i,j,iel,4) = shape*c_mat*c_mat ! p = rho * c^2 (acoustic)
102+
this%solution%interior(i,j,iel,5) = c_mat ! sound speed for this material
103+
enddo
104+
enddo
105+
enddo
106+
107+
call this%solution%UpdateDevice()
108+
109+
endsubroutine setInitialCondition
110+
111+
endmodule lineareuler2d_boneandmarrow_model
112+
113+
program LinearEuler_BoneAndMarrow
114+
115+
use self_data
116+
use lineareuler2d_boneandmarrow_model
117+
118+
implicit none
119+
120+
character(SELF_INTEGRATOR_LENGTH),parameter :: integrator = 'rk3'
121+
integer,parameter :: controlDegree = 3
122+
integer,parameter :: targetDegree = 6
123+
! CFL for DGSEM: dt < dx_min / (c_max * (N+1)^2). With curved elements
124+
! near the material interfaces dx_min can be much smaller than the
125+
! background 1.0, so we run conservatively.
126+
real(prec),parameter :: dt = 1.0e-3_prec
127+
real(prec),parameter :: endtime = 0.5_prec
128+
real(prec),parameter :: iointerval = 0.1_prec
129+
real(prec) :: e0,ef
130+
type(lineareuler2d_boneandmarrow) :: modelobj
131+
type(Lagrange),target :: interp
132+
type(Mesh2D),target :: mesh
133+
type(SEMQuad),target :: geometry
134+
character(LEN=255) :: WORKSPACE
135+
136+
! Read the multi-material ISM-MM mesh. Boundary names from the
137+
! .mesh file (e.g. "outer") become bc indices 1..nBCs; remap them
138+
! all to SELF_BC_RADIATION since the disk has one transparent
139+
! outer boundary.
140+
call get_environment_variable("WORKSPACE",WORKSPACE)
141+
call mesh%Read_HOHQMesh(trim(WORKSPACE)//"/share/mesh/MultiMaterial2D/BoneAndMarrow.mesh")
142+
call mesh%ResetBoundaryConditionType(SELF_BC_RADIATION)
143+
144+
call interp%Init(N=controlDegree, &
145+
controlNodeType=GAUSS, &
146+
M=targetDegree, &
147+
targetNodeType=UNIFORM)
148+
149+
call geometry%Init(interp,mesh%nElem)
150+
call geometry%GenerateFromMesh(mesh)
151+
152+
call modelobj%Init(mesh,geometry)
153+
modelobj%prescribed_bcs_enabled = .false.
154+
modelobj%tecplot_enabled = .false.
155+
modelobj%rho0 = 1.0_prec
156+
157+
call modelobj%setInitialCondition()
158+
159+
call modelobj%WriteModel()
160+
call modelobj%IncrementIOCounter()
161+
162+
call modelobj%CalculateEntropy()
163+
call modelobj%ReportEntropy()
164+
e0 = modelobj%entropy
165+
166+
call modelobj%SetTimeIntegrator(integrator)
167+
call modelobj%ForwardStep(endtime,dt,iointerval)
168+
ef = modelobj%entropy
169+
170+
if(ef /= ef) then
171+
print*,"Error: Final entropy is inf or nan",ef
172+
stop 1
173+
endif
174+
if(ef > 10.0_prec*e0) then
175+
print*,"Error: Final entropy grew unphysically (e0=",e0," ef=",ef,")"
176+
stop 1
177+
endif
178+
179+
call modelobj%free()
180+
call mesh%free()
181+
call geometry%free()
182+
call interp%free()
183+
184+
endprogram LinearEuler_BoneAndMarrow

‎examples/linear_euler2d_planewave_propagation.f90‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -47,6 +47,7 @@ module lineareuler2d_planewave_prop_model
4747
real(prec) :: x0 = 0.0_prec ! x component of the wave center position
4848
real(prec) :: y0 = 0.0_prec ! y component of the wave center position
4949
real(prec) :: L = 1.0_prec ! Halfwidth of the plane wave
50+
real(prec) :: c = 1.0_prec ! Reference sound speed (parameter for the analytical solution)
5051

5152
contains
5253

@@ -97,6 +98,7 @@ subroutine setInitialCondition(this)
9798
this%solution%interior(i,j,iel,2) = u*shape ! u
9899
this%solution%interior(i,j,iel,3) = v*shape ! v
99100
this%solution%interior(i,j,iel,4) = p*shape ! pressure
101+
this%solution%interior(i,j,iel,5) = this%c ! sound speed (uniform background)
100102

101103
enddo
102104

@@ -132,6 +134,7 @@ subroutine hbc2d_Prescribed_lineareuler2d_planewave(bc,mymodel)
132134
m%solution%extBoundary(i,j,iEl,2) = u*shape ! u
133135
m%solution%extBoundary(i,j,iEl,3) = v*shape ! v
134136
m%solution%extBoundary(i,j,iEl,4) = p*shape ! pressure
137+
m%solution%extBoundary(i,j,iEl,5) = m%c ! sound speed
135138
enddo
136139
enddo
137140
endselect

‎examples/linear_euler2d_planewave_reflection.f90‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -50,6 +50,7 @@ module lineareuler2d_planewave_model
5050
real(prec) :: x0 = 0.0_prec ! x component of the wave center position
5151
real(prec) :: y0 = 0.0_prec ! y component of the wave center position
5252
real(prec) :: L = 1.0_prec ! Halfwidth of the plane wave
53+
real(prec) :: c = 1.0_prec ! Reference sound speed (parameter for the analytical solution)
5354

5455
contains
5556

@@ -107,6 +108,7 @@ subroutine setInitialCondition(this)
107108
this%solution%interior(i,j,iel,2) = u*(shi-shr) ! u
108109
this%solution%interior(i,j,iel,3) = v*(shi+shr) ! v
109110
this%solution%interior(i,j,iel,4) = p*(shi+shr) ! pressure
111+
this%solution%interior(i,j,iel,5) = this%c ! sound speed (uniform background)
110112

111113
enddo
112114

@@ -147,6 +149,7 @@ subroutine hbc2d_Prescribed_lineareuler2d_planewave(bc,mymodel)
147149
m%solution%extBoundary(i,j,iEl,2) = u*(shi-shr) ! u
148150
m%solution%extBoundary(i,j,iEl,3) = v*(shi+shr) ! v
149151
m%solution%extBoundary(i,j,iEl,4) = p*(shi+shr) ! pressure
152+
m%solution%extBoundary(i,j,iEl,5) = m%c ! sound speed
150153
enddo
151154
enddo
152155
endselect

‎examples/linear_euler2d_spherical_soundwave_closeddomain.f90‎

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -73,10 +73,10 @@ program LinearEuler_Example
7373
modelobj%prescribed_bcs_enabled = .false. ! Disables prescribed boundary condition block for gpu accelerated implementations
7474
modelobj%tecplot_enabled = .false. ! Disables tecplot output
7575
modelobj%rho0 = rho0 ! optional, set the reference density
76-
modelobj%c = c ! optional set the reference sound wave speed
7776

78-
! Set the initial condition
79-
call modelobj%SphericalSoundWave(rhoprime,Lr,x0,y0)
77+
! Set the initial condition. The sound speed `c` is now carried as a
78+
! solution variable, so it is passed in here and written into solution(...,5).
79+
call modelobj%SphericalSoundWave(rhoprime,Lr,x0,y0,c)
8080

8181
call modelobj%WriteModel()
8282
call modelobj%IncrementIOCounter()

‎share/mesh/MultiMaterial2D/BoneAndMarrow.control‎

Lines changed: 97 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,97 @@
1+
%
2+
% ------------------------------------------------------------------
3+
% Control file for two materials embedded in a background material.
4+
% ------------------------------------------------------------------
5+
%
6+
% -------------
7+
% Control Block
8+
% -------------
9+
%
10+
\begin{CONTROL_INPUT}
11+
%
12+
% Plot and stats file names can be "none" if no output
13+
% is desired.
14+
%
15+
\begin{RUN_PARAMETERS}
16+
mesh file name = BoneAndMarrow.mesh
17+
plot file name = BoneAndMarrow.tec
18+
stats file name = BoneAndMarrow.txt
19+
mesh file format = ISM-MM
20+
polynomial order = 4
21+
plot file format = sem
22+
\end{RUN_PARAMETERS}
23+
24+
\begin{MESH_PARAMETERS}
25+
element type = quad
26+
\end{MESH_PARAMETERS}
27+
28+
\begin{BACKGROUND_GRID}
29+
background grid size = [1.0,1.0,0.0]
30+
\end{BACKGROUND_GRID}
31+
%
32+
% Spring smoother parameters are pretty standard after lots
33+
% of trials.
34+
%
35+
\begin{SPRING_SMOOTHER}
36+
smoothing type = LinearAndCrossBarSpring
37+
number of iterations = 25
38+
\end{SPRING_SMOOTHER}
39+
40+
\end{CONTROL_INPUT}
41+
%
42+
% -----------
43+
% Model Block
44+
% -----------
45+
%
46+
\begin{MODEL}
47+
48+
\begin{MATERIALS}
49+
background material name = Muscle
50+
\end{MATERIALS}
51+
%
52+
% If an outer boundary is not specified, then it is assumed to be a box
53+
% the outer boundary is implicitly a CHAIN. There is only one. Curves
54+
% implemented are SPLINE_CURVE, END_POINTS_LINE, PARAMETRIC_EQUATION_CURVE
55+
%
56+
\begin{OUTER_BOUNDARY}
57+
58+
\begin{PARAMETRIC_EQUATION_CURVE}
59+
name = outer
60+
xEqn = x(t) = 14.0*cos(2*pi*t)
61+
yEqn = y(t) = 14.0*sin(2*pi*t)
62+
zEqn = z(t) = 0.0
63+
\end{PARAMETRIC_EQUATION_CURVE}
64+
65+
\end{OUTER_BOUNDARY}
66+
67+
%
68+
% There are an arbitrary number of inner boundaries.
69+
% Each is a chain of curves, even if there is only
70+
% one curve in the chain.
71+
%
72+
\begin{INTERFACE_BOUNDARIES}
73+
74+
\begin{CHAIN}
75+
name = Bone
76+
\begin{PARAMETRIC_EQUATION_CURVE}
77+
name = Circle1
78+
xEqn = f(t) = 6.5*cos(2*pi*t) - 0.5
79+
yEqn = f(t) = 6.5*sin(2*pi*t)
80+
zEqn = z(t) = 0.0
81+
\end{PARAMETRIC_EQUATION_CURVE}
82+
\end{CHAIN}
83+
84+
\begin{CHAIN}
85+
name = Marrow
86+
\begin{PARAMETRIC_EQUATION_CURVE}
87+
name = Circle2
88+
xEqn = f(t) = -1.5 + 1.5*cos(2*pi*t)
89+
yEqn = f(t) = 1.5*sin(2*pi*t) + 0.5
90+
zEqn = z(t) = 0.0
91+
\end{PARAMETRIC_EQUATION_CURVE}
92+
\end{CHAIN}
93+
94+
\end{INTERFACE_BOUNDARIES}
95+
96+
\end{MODEL}
97+
\end{FILE}

0 commit comments

Comments
 (0)