magnetic sim agreees with calculations now
This commit is contained in:
jaseg
parent
a2d1429036
commit
cc59a6567e
4 changed files with 96 additions and 56 deletions
|
|
@ -9,31 +9,20 @@ ro4003c:
|
|||
Density: 1790 # 23°C
|
||||
Relative Permeability: 1
|
||||
Relative Permittivity: 3.55
|
||||
fr4:
|
||||
Density: 1850 # 23°C
|
||||
Relative Permeability: 1
|
||||
Relative Permittivity: 4.4
|
||||
Heat Conductivity: 0.81 # in-plane
|
||||
ideal:
|
||||
Relative Permittivity: 1
|
||||
copper:
|
||||
Density: 8960.0 # 0°C
|
||||
Electric Conductivity: 32300000 # 200°C
|
||||
Electric Conductivity: 59600000 # 20°C
|
||||
Emissivity: 0.012 # 327°C
|
||||
Heat Capacity: 415.0 # 200°C
|
||||
Heat Conductivity: 401.0 # 0°C
|
||||
Relative Permeability: 1
|
||||
Relative Permittivity: 1
|
||||
graphite_CZ3-R6300: # crucible
|
||||
Density: 1730.0
|
||||
Electric Conductivity: 58800
|
||||
Emissivity: 0.81 # 205°C
|
||||
Heat Capacity: 1237.0
|
||||
Heat Conductivity: 65 # 20°C
|
||||
Relative Permeability: 1
|
||||
Relative Permittivity: 1
|
||||
graphite_FU8957: # heater
|
||||
Density: 1750.0
|
||||
Emissivity: 0.81 # 250°C
|
||||
Heat Capacity: 1237.0
|
||||
Heat Conductivity: 105 # averaged over different given values
|
||||
Relative Permeability: 1
|
||||
Relative Permittivity: 1
|
||||
steel_1.4541:
|
||||
Density: 7900.0 # 20°C
|
||||
Electric Conductivity: 1370
|
||||
|
|
@ -42,26 +31,6 @@ steel_1.4541:
|
|||
Heat Conductivity: 15.0 # 20°C
|
||||
Relative Permeability: 1
|
||||
Relative Permittivity: 1
|
||||
tin_liquid:
|
||||
Density: 6980.0
|
||||
Electric Conductivity: 2080000
|
||||
Emissivity: 0.064 # set equal to solid
|
||||
Heat Capacity: 252.7
|
||||
Heat Conductivity: 29.0
|
||||
Relative Permeability: 1
|
||||
Relative Permittivity: 1
|
||||
Liquid: 'Logical True'
|
||||
tin_solid:
|
||||
Density: 7179.0
|
||||
Electric Conductivity: 4380000
|
||||
Emissivity: 0.064
|
||||
Heat Capacity: 244.0
|
||||
Heat Conductivity: 60.0
|
||||
Relative Permeability: 1
|
||||
Relative Permittivity: 1
|
||||
Solid: 'Logical True'
|
||||
Melting Point: 505
|
||||
Latent Heat: 59600
|
||||
water:
|
||||
Density: 1000.0
|
||||
Heat Capacity: 4182.0
|
||||
|
|
|
|||
|
|
@ -1,16 +1,17 @@
|
|||
Static_Current_Conduction:
|
||||
Equation: Static Current Conduction
|
||||
Variable: Potential
|
||||
Variable DOFs: 1
|
||||
Procedure: '"StatCurrentSolve" "StatCurrentSolver"'
|
||||
Calculate Volume Current: True
|
||||
Calculate Joule Heating: False
|
||||
Optimize Bandwidth: True
|
||||
Nonlinear System Max Iterations: 1
|
||||
Linear System Solver: Iterative
|
||||
Linear System Iterative Method: BiCGStabl
|
||||
Linear System Iterative Method: CG
|
||||
Linear System Max Iterations: 5000
|
||||
Linear System Convergence Tolerance: 1.0e-10
|
||||
Linear System Preconditioning: ILU0
|
||||
Linear System Preconditioning: ILU3
|
||||
Linear System ILUT Tolerance: 1.0e-3
|
||||
Linear System Abort Not Converged: False
|
||||
Linear System Residual Output: 20
|
||||
|
|
@ -21,6 +22,7 @@ Magneto_Dynamics:
|
|||
Procedure: '"MagnetoDynamics" "WhitneyAVSolver"'
|
||||
Newton-Raphson Iteration: True
|
||||
Nonlinear System Max Iterations: 1
|
||||
Fix Input Current Density: True
|
||||
Linear System Solver: Iterative
|
||||
Linear System Preconditioning: none
|
||||
Linear System Convergence Tolerance: 1e-8
|
||||
|
|
|
|||
|
|
@ -1,15 +1,16 @@
|
|||
#!/usr/bin/env python3
|
||||
|
||||
import math
|
||||
from pathlib import Path
|
||||
import multiprocessing
|
||||
import re
|
||||
import tempfile
|
||||
import subprocess
|
||||
import fnmatch
|
||||
import shutil
|
||||
import numpy as np
|
||||
|
||||
from pyelmer import elmer
|
||||
import subprocess_tee
|
||||
import click
|
||||
from scipy import constants
|
||||
|
||||
|
|
@ -31,6 +32,17 @@ def enumerate_mesh_bodies(msh_file):
|
|||
tag, _, name = line.partition(' ')
|
||||
yield name.strip().strip('"'), (int(dim), int(tag))
|
||||
|
||||
# https://en.wikipedia.org/wiki/Metric_prefix
|
||||
SI_PREFIX = 'QRYZEPTGMk mµnpfazyrq'
|
||||
|
||||
def format_si(value, unit='', fractional_digits=1):
|
||||
mag = int(math.log10(abs(value))//3)
|
||||
value /= 1000**mag
|
||||
prefix = SI_PREFIX[SI_PREFIX.find(' ') - mag].strip()
|
||||
value = f'{{:.{fractional_digits}f}}'.format(value)
|
||||
return f'{value} {prefix}{unit}'
|
||||
|
||||
|
||||
INPUT_EXT_MAP = {
|
||||
'.grd': 1,
|
||||
'.mesh*': 2,
|
||||
|
|
@ -56,7 +68,7 @@ OUTPUT_EXT_MAP = {
|
|||
'.msh': 4,
|
||||
'.vtu': 5}
|
||||
|
||||
def elmer_grid(infile, outfile=None, intype=None, outtype=None, cwd=None, **kwargs):
|
||||
def elmer_grid(infile, outfile=None, intype=None, outtype=None, cwd=None, stdout_log=None, stderr_log=None, **kwargs):
|
||||
infile = Path(infile)
|
||||
if outfile is not None:
|
||||
outfile = Path(outfile)
|
||||
|
|
@ -80,10 +92,20 @@ def elmer_grid(infile, outfile=None, intype=None, outtype=None, cwd=None, **kwar
|
|||
args.extend(str(v) for v in value)
|
||||
else:
|
||||
args.append(str(value))
|
||||
subprocess.run(args, cwd=cwd, check=True)
|
||||
|
||||
def elmer_solver(cwd):
|
||||
subprocess.run(['ElmerSolver'], cwd=cwd)
|
||||
result = subprocess_tee.run(args, cwd=cwd, check=True)
|
||||
if stdout_log:
|
||||
Path(stdout_log).write_text(result.stdout or '')
|
||||
if stderr_log:
|
||||
Path(stderr_log).write_text(result.stderr or '')
|
||||
|
||||
def elmer_solver(cwd, stdout_log=None, stderr_log=None):
|
||||
result = subprocess_tee.run(['ElmerSolver'], cwd=cwd, check=True)
|
||||
if stdout_log:
|
||||
Path(stdout_log).write_text(result.stdout or '')
|
||||
if stderr_log:
|
||||
Path(stderr_log).write_text(result.stderr or '')
|
||||
return result
|
||||
|
||||
|
||||
@click.command()
|
||||
|
|
@ -132,14 +154,17 @@ def run_capacitance_simulation(mesh_file, sim_dir):
|
|||
boundary_airbox.data['Electric Infinity BC'] = 'True'
|
||||
|
||||
with tempfile.TemporaryDirectory() as tmpdir:
|
||||
if sim_dir:
|
||||
tmpdir = str(sim_dir)
|
||||
tmpdir = sim_dir if sim_dir else Path(tmpdir)
|
||||
|
||||
sim.write_startinfo(tmpdir)
|
||||
sim.write_sif(tmpdir)
|
||||
# Convert mesh from gmsh to elemer formats. Also scale it from 1 unit = 1 mm to 1 unit = 1 m (SI units)
|
||||
elmer_grid(mesh_file.absolute(), 'mesh', cwd=tmpdir, scale=[1e-3, 1e-3, 1e-3])
|
||||
elmer_solver(tmpdir)
|
||||
elmer_grid(mesh_file.absolute(), 'mesh', cwd=tmpdir, scale=[1e-3, 1e-3, 1e-3],
|
||||
stdout_log=(tmpdir / 'ElmerGrid_stdout.log'),
|
||||
stderr_log=(tmpdir / 'ElmerGrid_stderr.log'))
|
||||
elmer_solver(tmpdir,
|
||||
stdout_log=(tmpdir / 'ElmerSolver_stdout.log'),
|
||||
stderr_log=(tmpdir / 'ElmerSolver_stderr.log'))
|
||||
|
||||
capacitance_matrix = np.loadtxt(tmpdir / 'capacitance.txt')
|
||||
|
||||
|
|
@ -194,6 +219,9 @@ def run_inductance_simulation(mesh_file, sim_dir):
|
|||
bdy_if_bottom.material = copper
|
||||
bdy_if_bottom.equation = copper_eqn
|
||||
|
||||
potential_force = elmer.BodyForce(sim, 'electric_potential', {'Electric Potential': 'Equals "Potential"'})
|
||||
bdy_trace.body_force = potential_force
|
||||
|
||||
# boundaries
|
||||
boundary_airbox = elmer.Boundary(sim, 'FarField', [physical['airbox_surface'][1]])
|
||||
boundary_airbox.data['Electric Infinity BC'] = 'True'
|
||||
|
|
@ -206,15 +234,46 @@ def run_inductance_simulation(mesh_file, sim_dir):
|
|||
boundary_vminus.data['Potential'] = 0.0
|
||||
|
||||
with tempfile.TemporaryDirectory() as tmpdir:
|
||||
if sim_dir:
|
||||
tmpdir = str(sim_dir)
|
||||
tmpdir = sim_dir if sim_dir else Path(tmpdir)
|
||||
|
||||
sim.write_startinfo(tmpdir)
|
||||
sim.write_sif(tmpdir)
|
||||
# Convert mesh from gmsh to elemer formats. Also scale it from 1 unit = 1 mm to 1 unit = 1 m (SI units)
|
||||
elmer_grid(mesh_file.absolute(), 'mesh', cwd=tmpdir, scale=[1e-3, 1e-3, 1e-3])
|
||||
elmer_solver(tmpdir)
|
||||
elmer_grid(mesh_file.absolute(), 'mesh', cwd=tmpdir, scale=[1e-3, 1e-3, 1e-3],
|
||||
stdout_log=(tmpdir / 'ElmerGrid_stdout.log'),
|
||||
stderr_log=(tmpdir / 'ElmerGrid_stderr.log'))
|
||||
solver_stdout, solver_stderr = (tmpdir / 'ElmerSolver_stdout.log'), (tmpdir / 'ElmerSolver_stderr.log')
|
||||
res = elmer_solver(tmpdir,
|
||||
stdout_log=solver_stdout,
|
||||
stderr_log=solver_stderr)
|
||||
|
||||
P, R, U_mag = None, None, None
|
||||
solver_error = False
|
||||
for l in res.stdout.splitlines():
|
||||
if (m := re.fullmatch(r'StatCurrentSolve:\s*Total Heating Power\s*:\s*([0-9.+-Ee]+)\s*', l)):
|
||||
P = float(m.group(1))
|
||||
elif (m := re.fullmatch(r'StatCurrentSolve:\s*Effective Resistance\s*:\s*([0-9.+-Ee]+)\s*', l)):
|
||||
R = float(m.group(1))
|
||||
elif (m := re.fullmatch(r'MagnetoDynamicsCalcFields:\s*ElectroMagnetic Field Energy\s*:\s*([0-9.+-Ee]+)\s*', l)):
|
||||
U_mag = float(m.group(1))
|
||||
elif re.fullmatch(r'IterSolve: Linear iteration did not converge to tolerance', l):
|
||||
solver_error = True
|
||||
|
||||
if solver_error:
|
||||
raise click.ClickException(f'Error: One of the solvers did not converge. See log files for details:\n{solver_stdout.absolute()}\n{solver_stderr.absolute()}')
|
||||
elif P is None or R is None or U_mag is None:
|
||||
raise click.ClickException(f'Error during solver execution. Electrical parameters could not be calculated. See log files for details:\n{solver_stdout.absolute()}\n{solver_stderr.absolute()}')
|
||||
|
||||
U = math.sqrt(P*R)
|
||||
I = math.sqrt(P/R)
|
||||
L = 2*U_mag / (I**2)
|
||||
|
||||
assert math.isclose(U, 1.0, abs_tol=1e-3)
|
||||
|
||||
print(f'Coil resistance calculated by solver: {format_si(R, "Ω")}')
|
||||
print(f'Inductance calucated from field: {format_si(L, "H")}')
|
||||
|
||||
|
||||
if __name__ == '__main__':
|
||||
run_inductance_simulation()
|
||||
|
||||
|
|
|
|||
|
|
@ -317,6 +317,8 @@ def traces_to_gmsh_mag(traces, mesh_out, bbox, model_name='gerbonara_board', log
|
|||
airbox_physical = gmsh.model.add_physical_group(3, [airbox], name='airbox')
|
||||
trace_physical = gmsh.model.add_physical_group(3, [toplevel_tag], name='trace')
|
||||
|
||||
gmsh.model.mesh.setSize([(0, tag) for dim, tag in gmsh.model.getBoundary([(3, toplevel_tag)], recursive=True) if dim == 0], 0.300)
|
||||
|
||||
#interface_tags_top = gmsh.model.getBoundary([(3, contact_tag_top)], oriented=False)
|
||||
#interface_tags_bottom = gmsh.model.getBoundary([(3, contact_tag_bottom)], oriented=False)
|
||||
|
||||
|
|
@ -464,6 +466,8 @@ def print_valid_twists(ctx, param, value):
|
|||
@click.option('--via-offset', type=float, default=None, help='Radially offset vias from trace endpoints [mm]')
|
||||
@click.option('--keepout-zone/--no-keepout-zone', default=True, help='Add a keepout are to the footprint (default: yes)')
|
||||
@click.option('--keepout-margin', type=float, default=5, help='Margin between outside of coil and keepout area (mm, default: 5)')
|
||||
@click.option('--copper-thickness', type=float, default=0.035, help='Copper thickness for resistance calculation and mesh generation in mm. Default: 0.035mm ^= 1 Oz')
|
||||
@click.option('--board-thickness', type=float, default=1.53, help='Board substrate thickness for mesh generation in mm. Default: 1.53mm')
|
||||
@click.option('--twists', type=int, default=1, help='Number of twists per revolution. Note that this number must be co-prime to the number of turns. Run with --show-twists to list valid values. (default: 1)')
|
||||
@click.option('--circle-segments', type=int, default=64, help='When not using arcs, the number of points to use for arc interpolation per 360 degrees.')
|
||||
@click.option('--show-twists', callback=print_valid_twists, expose_value=False, type=int, is_eager=True, help='Calculate and show valid --twists counts for the given number of turns. Takes the number of turns as a value.')
|
||||
|
|
@ -477,7 +481,8 @@ def print_valid_twists(ctx, param, value):
|
|||
@click.version_option()
|
||||
def generate(outfile, turns, outer_diameter, inner_diameter, via_diameter, via_drill, via_offset, trace_width, clearance,
|
||||
footprint_name, layer_pair, twists, clipboard, counter_clockwise, keepout_zone, keepout_margin,
|
||||
arc_tolerance, pcb, mesh_out, magneticalc_out, circle_segments, mag_mesh_out):
|
||||
arc_tolerance, pcb, mesh_out, magneticalc_out, circle_segments, mag_mesh_out, copper_thickness,
|
||||
board_thickness):
|
||||
if 'WAYLAND_DISPLAY' in os.environ:
|
||||
copy, paste, cliputil = ['wl-copy'], ['wl-paste'], 'xclip'
|
||||
else:
|
||||
|
|
@ -722,7 +727,12 @@ def generate(outfile, turns, outer_diameter, inner_diameter, via_diameter, via_d
|
|||
svg_vias.append(Tag('circle', cx=xv, cy=yv, r=via_diameter/2, stroke='none', fill='white'))
|
||||
svg_vias.append(Tag('circle', cx=xv, cy=yv, r=via_drill/2, stroke='none', fill='black'))
|
||||
|
||||
print(f'Approximate track length: {clen*twists*2:.2f} mm', file=sys.stderr)
|
||||
l_total = clen*twists*2
|
||||
print(f'Approximate track length: {l_total:.2f} mm', file=sys.stderr)
|
||||
A = copper_thickness/1e3 * trace_width/1e3
|
||||
rho = 1.68e-8
|
||||
R = l_total/1e3 * rho / A
|
||||
print(f'Approximate resistance: {R:g} Ω', file=sys.stderr)
|
||||
|
||||
top_pad = make_pad(1, [layer_pair[0]], outer_radius, 0)
|
||||
pads.append(top_pad)
|
||||
|
|
@ -812,10 +822,10 @@ def generate(outfile, turns, outer_diameter, inner_diameter, via_diameter, via_d
|
|||
|
||||
r = outer_diameter/2 + 20
|
||||
if mesh_out:
|
||||
traces_to_gmsh(traces, mesh_out, ((-r, -r), (r, r)))
|
||||
traces_to_gmsh(traces, mesh_out, ((-r, -r), (r, r)), copper_thickness=copper_thickness, board_thickness=board_thickness)
|
||||
|
||||
if mag_mesh_out:
|
||||
traces_to_gmsh_mag(traces, mag_mesh_out, ((-r, -r), (r, r)))
|
||||
traces_to_gmsh_mag(traces, mag_mesh_out, ((-r, -r), (r, r)), copper_thickness=copper_thickness, board_thickness=board_thickness)
|
||||
|
||||
if magneticalc_out:
|
||||
traces_to_magneticalc(traces, magneticalc_out)
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue