Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Use spline surface fitting for in-vessel component CAD generation #171

Draft
wants to merge 6 commits into
base: main
Choose a base branch
from

Conversation

connoramoreno
Copy link
Collaborator

@connoramoreno connoramoreno commented Oct 25, 2024

Overview

Modifies the in-vessel component CAD generation workflow to use CadQuery’s spline surface fitting routine, instead of its lofting routine, to define the outer surfaces of each in-vessel component (IVC). The benefits of this change are several-fold:

  • Facilitates radial build definitions that begin at toroidal and poloidal angles other than 0 degrees
  • Builds can span any toroidal extent – no longer need to repeat periods
  • The spline surface fitting routine tends to be more robust against steep gradients in IVC thickness matrices
  • IVC outer surfaces tend to be smoother
  • IVC file sizes are reduced

List of changes in branch:

  • Modifies ParaStell’s Surface.generate_surface method to use CadQuery’s Face.makeSplineApprox surface spline technique instead of its Solid.makeLoft routine
  • Modifies ParaStell’s InVesselBuild.generate_components method to construct CAD solids by assembling necessary surfaces, instead of performing a cut operation
  • Adds InVesselBuild._check_edge method
  • Removes the now extraneous Rib class and moves any necessary methods to the Surface class
  • Renames num_ribs and num_rib_pts parameters due to removal of Rib class
  • Removes the now extraneous repeat parameter
  • Introduces new parameter transpose_fit, which is a flag to indicate that the 2-D iterable of points through which each spline surface is fit should be transposed. This swaps the order in which the grid axes are iterated over. For a yet-to-be-understood reason, this can fix errant CAD generation.
  • Adds min_mesh_size and max_mesh_size arguments to the InVesselBuild.export_cad_to_dagmc method
  • Modifies the utils.normalize function to allow 3-D NumPy arrays
  • Updates documentation

Issues to be solved:

  • While CAD generation more robust than lofting workflow, still temperamental to some degree
  • In-vessel build test failing for unusual reason

Closes #26
Closes #170

Illustration of new CAD generation workflow:

Screenshot 2024-10-25 at 11 53 47 AM

Figure 1. Create separate surfaces for outer surface and end caps (plasma/chamber).

Screenshot 2024-10-25 at 11 54 01 AM

Figure 2. Combine surfaces into single solid (plasma/chamber).

Screenshot 2024-10-25 at 11 54 32 AM

Figure 3. Create separate surfaces for inner surface (outer surface of interior component), outer surface, and end caps.

Screenshot 2024-10-25 at 11 54 46 AM

Figure 4. Combine surfaces into single solid.

Demonstration of new capabilities:

Screenshot 2024-10-25 at 12 34 27 PM

Figure 5. 360-degree geometry.

Screenshot 2024-10-25 at 12 40 15 PM

Figure 6. 360-degree geometry (cut at mid-plane, visualized in ParaView).

Screenshot 2024-10-25 at 11 42 07 AM

Figure 7. In-vessel build with thickness matrix of component defined as 200 * numpy.random.rand(9, 9).

@Edgar-21
Copy link
Contributor

I'd just like to mention that I started using this this weekend and have found it to be very useful, thanks @connoramoreno

Copy link
Member

@gonuke gonuke left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A few minor changes, I think. Do we need to update examples, as well?

Will the change in interface break all the old parastell usage?

Comment on lines +236 to +237
"""if self.radial_build.toroidal_angles[0] >= 0.0:
toroidal_angle = (toroidal_angle + 2 * np.pi) % (2 * np.pi)"""
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we need this? Seems like we might want it, in general?

@@ -226,39 +217,71 @@ def calculate_loci(self):

[surface.calculate_loci() for surface in self.Surfaces.values()]

def _check_edge(self, edge):
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Functions returning booleans may be better with a name that indicates the outcome, like _is_near_end

Comment on lines +260 to +265
edge_list = [
edge
for edge in outer_surface.Edges()
if self._check_edge(edge)
]
wire_list = [cq.Wire.assembleEdges([edge]) for edge in edge_list]
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can't we just combine these? I don't see edge_list used elsewhere:

Suggested change
edge_list = [
edge
for edge in outer_surface.Edges()
if self._check_edge(edge)
]
wire_list = [cq.Wire.assembleEdges([edge]) for edge in edge_list]
wire_list = [cq.Wire.assembleEdges([edge]) for edge in outer_surface.Edges()
if self._check_edge(edge)]

@@ -460,8 +446,11 @@ def _vmec2xyz(self, poloidal_offset=0):
"""
return self.scale * np.array(
[
self.vmec_obj.vmec2xyz(self.s, theta, self.phi)
for theta in (self.theta_list + poloidal_offset)
[
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Need to update the docstring to match this new behavior

@@ -619,14 +620,15 @@ def poloidal_angles(self, angle_list):
self._logger.error(e.args[0])
raise e

self._poloidal_angles = angle_list
if self._poloidal_angles[-1] - self._poloidal_angles[0] > 360.0:
if angle_list[-1] - angle_list[0] != 360.0:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Exactly equal in floating point math is tricky.

Should we check that these are close and then move one to make them exact if they are close, or make an error if they aren't close?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
3 participants