Multiphysics: HFSS-Mechanical MRI analysis

The goal of this workshop is to use a coil tuned to 63.8 MHz to determine the temperature rise in a gel phantom near an implant given a background SAR of 1 W/kg.

Steps to follow Step 1: Simulate coil loaded by empty phantom: Scale input to coil ports to produce desired background SAR of 1 W/kg at location that will later contain the implant. Step 2: Simulate coil loaded by phantom containing implant in proper location: View SAR in tissue surrounding implant. Step 3: Thermal simulation: Link HFSS to transient thermal solver to find temperature rise in tissue near implant vs. time.

Perform required imports

Perform required imports.

import os.path

from pyaedt import Hfss, Mechanical, Icepak, downloads

Set AEDT version

Set AEDT version.

aedt_version = "2024.1"

Set non-graphical mode

Set non-graphical mode. ` You can set non_graphical either to True or False.

non_graphical = False

Project load

Open the ANSYS Electronics Desktop 2018.2 Open project background_SAR.aedt Project contains phantom and airbox Phantom consists of two objects: phantom and implant_box Separate objects are used to selectively assign mesh operations Material properties defined in this project already contain #electrical and thermal properties.

project_path = downloads.download_file(directory="mri")
hfss = Hfss(os.path.join(project_path, "background_SAR.aedt"), specified_version=aedt_version, non_graphical=non_graphical,
            new_desktop_session=True)

C:\actions-runner\_work\_tool\Python\3.10.9\x64\lib\subprocess.py:1072: ResourceWarning: subprocess 8144 is still running
  _warn("subprocess %s is still running" % self.pid,
C:\actions-runner\_work\pyaedt\pyaedt\testenv\lib\site-packages\pyaedt\generic\settings.py:383: ResourceWarning: unclosed file <_io.TextIOWrapper name='D:\\Temp\\pyaedt_ansys.log' mode='a' encoding='cp1252'>
  self._logger = val

Insert 3D component

The MRI Coil is saved as a separate 3D Component ‒ 3D Components store geometry (including parameters), material properties, boundary conditions, mesh assignments, and excitations ‒ 3D Components make it easy to reuse and share parts of a simulation

hfss.modeler.insert_3d_component(os.path.join(project_path, "coil.a3dcomp"))

<pyaedt.modeler.cad.components_3d.UserDefinedComponent object at 0x00000279A4C8D810>

Expression Cache

On the expression cache tab, define additional convergence criteria for self impedance of the four coil

ports ‒ Set each of these convergence criteria to 2.5 ohm For this demo number of passes is limited to 2 to reduce simulation time.

im_traces = hfss.get_traces_for_plot(get_mutual_terms=False, category="im(Z", first_element_filter="Coil1_p*")

hfss.setups[0].enable_expression_cache(
    report_type="Modal Solution Data",
    expressions=im_traces,
    isconvergence=True,
    isrelativeconvergence=False,
    conv_criteria=2.5,
    use_cache_for_freq=False)
hfss.setups[0].props["MaximumPasses"] = 2

Edit Sources

The 3D Component of the MRI Coil contains all the ports, but the sources for these ports are not yet defined. Browse to and select sources.csv. These sources were determined by tuning this coil at 63.8 MHz. Notice the “*input_scale” multiplier to allow quick adjustment of the coil excitation power.

hfss.edit_sources_from_file(os.path.join(project_path, "sources.csv"))

True

Run Simulation

Save and analyze the project.

hfss.save_project(os.path.join(project_path, "solved.aedt"))
hfss.analyze(cores=6)

True

Plot SAR on cut plane in phantom

Ensure that the SAR averaging method is set to Gridless Plot averagedSAR on GlobalYZ plane Draw Point1 at origin of the implant coordinate system

hfss.sar_setup(-1, tissue_mass=1, material_density=1, average_sar_method=1)
hfss.post.create_fieldplot_cutplane(assignment="implant:YZ", quantity="Average_SAR", filter_objects=["implant_box"])

hfss.modeler.set_working_coordinate_system("implant")
hfss.modeler.create_point([0, 0, 0], name="Point1")

hfss.post.plot_field(quantity="Average_SAR", assignment="implant:YZ", plot_type="CutPlane", show=False,
                     show_legend=False, filter_objects=["implant_box"])

<pyaedt.generic.plot.ModelPlotter object at 0x00000279A5899870>

Adjust Input Power to MRI Coil

The goal is to adjust the MRI coil’s input power, so that the averageSAR at Point1 is 1 W/kg Note that SAR and input power are linearly related To determine required input, calculate input_scale = 1/AverageSAR at Point1

sol_data = hfss.post.get_solution_data(expressions="Average_SAR",
                                       primary_sweep_variable="Freq",
                                       context="Point1",
                                       report_category="Fields")
sol_data.data_real()

hfss["input_scale"] = 1 / sol_data.data_real()[0]

Phantom with Implant

Import implant geometry. Subtract rod from implant_box. Assign titanium to the imported object rod. Analyze the project.

hfss.modeler.import_3d_cad(os.path.join(project_path, "implant_rod.sat"))

hfss.modeler["implant_box"].subtract("rod", keep_originals=True)
hfss.modeler["rod"].material_name = "titanium"
hfss.analyze(cores=6)
hfss.save_project()

True

Thermal Simulation

Initialize a new Mechanical Transient Thermal analysis. Mechanical Transient Thermal is available in AEDT from 2023 R2 as a Beta feature.

mech = Mechanical(solution_type="Transient Thermal", specified_version=aedt_version)

Copy geometries

Copy bodies from the HFSS project. 3D Component will not be copied.

mech.copy_solid_bodies_from(hfss)

True

Link sources to EM losses

Link sources to the EM losses. Assign external convection.

exc = mech.assign_em_losses(design=hfss.design_name, setup=hfss.setups[0].name, sweep="LastAdaptive",
                            map_frequency=hfss.setups[0].props["Frequency"],
                            surface_objects=mech.get_all_conductors_names())
mech.assign_uniform_convection(mech.modeler["Region"].faces, convection_value=1)

<pyaedt.modules.Boundary.BoundaryObject object at 0x00000279A1D93BE0>

Create Setup

Create a new setup and edit properties. Simulation will be for 60 seconds.

setup = mech.create_setup()
# setup.add_mesh_link("backgroundSAR")
# mech.create_dataset1d_design("PowerMap", [0, 239, 240, 360], [1, 1, 0, 0])
# exc.props["LossMultiplier"] = "pwl(PowerMap,Time)"

mech.modeler.set_working_coordinate_system("implant")
mech.modeler.create_point([0, 0, 0], name="Point1")
setup.props["Stop Time"] = 60
setup.props["Time Step"] = "10s"
setup.props["SaveFieldsType"] = "Every N Steps"
setup.props["N Steps"] = "2"

Analyze Mechanical

Analyze the project.

mech.analyze(cores=6)

True

Plot fields

Plot Temperature on cut plane. Plot Temperature on point.

mech.post.create_fieldplot_cutplane("implant:YZ", "Temperature", filter_objects=["implant_box"])
mech.save_project()

data = mech.post.get_solution_data("Temperature", primary_sweep_variable="Time", context="Point1",
                                   report_category="Fields")
#data.plot()

mech.post.plot_animated_field(quantity="Temperature", assignment="implant:YZ", plot_type="CutPlane",
                              intrinsics={"Time": "10s"}, variation_variable="Time",
                              variations=["10s", "20s", "30s", "40s", "50s", "60s"],
                              show=False, filter_objects=["implant_box"])

<pyaedt.generic.plot.ModelPlotter object at 0x00000279A0B1A320>

Thermal Simulation

Initialize a new Icepak Transient Thermal analysis.

ipk = Icepak(solution_type="Transient", specified_version=aedt_version)
ipk.design_solutions.problem_type = "TemperatureOnly"

Copy geometries

Copy bodies from the HFSS project. 3D Component will not be copied.

ipk.modeler.delete("Region")
ipk.copy_solid_bodies_from(hfss)

True

Link sources to EM losses

Link sources to the EM losses. Assign external convection.

exc = ipk.assign_em_losses(design=hfss.design_name, setup=hfss.setups[0].name, sweep="LastAdaptive",
                           map_frequency=hfss.setups[0].props["Frequency"],
                           surface_objects=ipk.get_all_conductors_names())

Create Setup

Create a new setup and edit properties. Simulation will be for 60 seconds.

setup = ipk.create_setup()

setup.props["Stop Time"] = 60
setup.props["N Steps"] = 2
setup.props["Time Step"] = 5
setup.props['Convergence Criteria - Energy'] = 1e-12

Mesh Region

Create a new mesh region and change accuracy level to 4.

bound = ipk.modeler["implant_box"].bounding_box
mesh_box = ipk.modeler.create_box(bound[:3], [bound[3] - bound[0], bound[4] - bound[1], bound[5] - bound[2]])
mesh_box.model = False
mesh_region = ipk.mesh.assign_mesh_region([mesh_box.name])
mesh_region.UserSpecifiedSettings = False
mesh_region.Level = 4
mesh_region.update()

True

Point Monitor

Create a new point monitor.

ipk.modeler.set_working_coordinate_system("implant")
ipk.monitor.assign_point_monitor([0, 0, 0], monitor_name="Point1")
ipk.assign_openings(ipk.modeler["Region"].top_face_z)

<pyaedt.modules.Boundary.BoundaryObject object at 0x00000279A5B1AF50>

Analyze and plot fields

Analyze the project. Plot temperature on cut plane. Plot temperature on monitor point.

ipk.analyze(cores=4, tasks=4)
ipk.post.create_fieldplot_cutplane("implant:YZ", "Temperature", filter_objects=["implant_box"])
ipk.save_project()

data = ipk.post.get_solution_data("Point1.Temperature", primary_sweep_variable="Time", report_category="Monitor")
#data.plot()

ipk.release_desktop(True, True)

True