bifilar-electromagnet-1.irmf

/*{
  "author": "Glenn M. Lewis",
  "copyright": "Apache-2.0",
  "date": "2020-02-13",
  "irmf": "1.0",
  "materials": ["copper","copper","dielectric"],
  "max": [25,25,71],
  "min": [-25,-25,-63],
  "notes": "The IRMF shader that started it all... the bifilar electromagnet.",
  "options": {
    "resolution": 512,
    "color1": [255,0,0,1],
    "color2": [0,255,0,1],
    "color3": [0,0,255,1]
  },
  "title": "bifilar electromagnet",
  "units": "mm",
  "version": "1.0"
}*/

#include "github.com/gmlewis/irmf-examples/blob/master/examples/012-bifilar-electromagnet/rotation.glsl"
#include "github.com/gmlewis/irmf-examples/blob/master/examples/012-bifilar-electromagnet/primitives.glsl"

#define M_PI 3.1415926535897932384626433832795

float coilSquareFace(float radius, float size, float gap, float nTurns, float trimEndAngle, in vec3 xyz) {
  // First, trivial reject on the two vertical ends of the coil.
  if (xyz.z < -0.5 * size || xyz.z > nTurns * (size + gap) + 0.5 * size) { return 0.0; }
  
  // Then, constrain the coil to the cylinder with wall thickness "size":
  float rxy = length(xyz.xy);
  if (rxy < radius - 0.5 * size || rxy > radius + 0.5 * size) { return 0.0; }
  
  // If the current point is between the coils, return no material:
  float angle = atan(xyz.y, xyz.x) / (2.0 * M_PI);
  if (angle < 0.0) { angle += 1.0; } // 0 <= angle <= 1 between coils from center to center.
  // 0 <= dz <= (size+gap) between coils from center to center.
  float dz = mod(xyz.z, size + gap);
  
  float lastHelixZ = angle * (size + gap);
  float coilNum = 0.0;
  if (lastHelixZ > dz) {
    lastHelixZ -= (size + gap); // center of current coil.
    coilNum = -1.0;
  }
  float nextHelixZ = lastHelixZ + (size + gap); // center of next higher vertical coil.
  
  // If the current point is within the gap between the two coils, reject it.
  if (dz > lastHelixZ + 0.5 * size && dz < nextHelixZ - 0.5 * size) { return 0.0; }
  
  coilNum += floor((xyz.z + (0.5 * size) - lastHelixZ) / (size + gap));
  
  // If the current point is in a coil numbered outside the current range, reject it.
  if (coilNum < 0.0 || coilNum >= nTurns) { return 0.0; }
  
  // TODO: Incorporate the trimEndAngle.
  
  return 1.0;
}

float coilPlusConnectorWires(int coilNum, int numCoils, float inc, float innerRadius, float connectorRadius, float size, float gap, float nTurns, in vec3 xyz) {
  float radiusOffset = float(coilNum - 1);
  mat4 xfm = mat4(1) * rotZ(radiusOffset * inc);
  float coilRadius = radiusOffset + innerRadius;
  float trimEndAngle = 2.0 * inc;
  if (coilNum == numCoils) {
    trimEndAngle = 0.5 * inc; // Special case to access the exit wire.
  } else if (coilNum == numCoils - 1) {
    trimEndAngle = 3.0 * inc;
  }
  xyz = (vec4(xyz, 1.0) * xfm).xyz;
  float coil = coilSquareFace(coilRadius, size, gap, nTurns, trimEndAngle, xyz);
  
  float bz = -(size + gap);
  float tz = nTurns * (size + gap);
  float tzp1 = (nTurns + 1.0) * (size + gap);
  
  coil += box(vec3(coilRadius, 0.0, 0.0), vec3(coilRadius, 0.0, bz), size, xyz);
  coil += box(vec3(coilRadius, 0.0, bz), vec3(connectorRadius, 0.0, bz), size, xyz);
  coil += box(vec3(connectorRadius, 0.0, bz), vec3(connectorRadius, 0.0, tzp1), size, xyz);
  
  float zexit = (nTurns + 10.0) * (size + gap);
  if (coilNum >= 3) { // Connect the start of this coil to the end of two coils prior.
    float lastCoilRadius = radiusOffset - 2.0 + innerRadius;
    coil += box(vec3(lastCoilRadius, 0.0, tzp1), vec3(connectorRadius, 0.0, tzp1), size, xyz);
    coil += box(vec3(lastCoilRadius, 0.0, tz), vec3(lastCoilRadius, 0.0, tzp1), size, xyz);
  } else if (coilNum == 2) { // Connect the start of 2 to the end of the last odd coil.
    float endOddRadius = float(numCoils - 2) + innerRadius;
    coil += box(vec3(endOddRadius, 0.0, tzp1), vec3(connectorRadius, 0.0, tzp1), size, xyz);
    coil += box(vec3(endOddRadius, 0.0, tz), vec3(endOddRadius, 0.0, tzp1), size, xyz);
  } else if (coilNum == 1) { // Bring out the exit wires.
    // Start of coil1:
    coil += box(vec3(connectorRadius, 0.0, tz), vec3(connectorRadius, 0.0, zexit), size, xyz);
  }
  
  if (coilNum == numCoils) { // Special case to access the exit wire.
    // End of coil 'numCoils':
    xfm = mat4(1) * rotZ(0.5 * inc);
    xyz = (vec4(xyz, 1.0) * xfm).xyz;
    coil += box(vec3(connectorRadius - (size + gap), 0.0, tz), vec3(connectorRadius - (size + gap), 0.0, zexit), size, xyz);
  }
  
  return coil;
}

vec3 bifilarElectromagnet(int numPairs, float innerRadius, float size, float gap, int numTurns, in vec3 xyz) {
  float nTurns = float(numTurns);
  int numCoils = 2*numPairs;
  float inc = 2.0 * M_PI / float(numCoils);
  float connectorRadius = innerRadius + float(numCoils) * (size + gap);
  
  float metal1 = 0.0;
  float metal2 = 0.0;
  
  for(int i = 1; i <= numCoils; i += 2) {
    metal1 += coilPlusConnectorWires(i, numCoils, inc, innerRadius, connectorRadius, size, gap, nTurns, xyz);
    metal2 += coilPlusConnectorWires(i + 1, numCoils, inc, innerRadius, connectorRadius, size, gap, nTurns, xyz);
  }
  metal1 = clamp(metal1, 0.0, 1.0);
  metal2 = clamp(metal2, 0.0, 1.0);
  
  float dielectric = 0.0;
  float dielectricRadius = 2.0 * float(numPairs) * (size + gap) + innerRadius + gap;
  float dielectricHeight = (nTurns + 4.0) * (size + gap);
  dielectric += cylinder(dielectricRadius, dielectricHeight, xyz + vec3(0, 0, 2));
  
  float spindleRadius = 2.0 * float(numPairs + 1) * (size + gap) + innerRadius;
  dielectric += cylinder(spindleRadius, 2.0 * (size + gap), xyz + vec3(0, 0, 2));
  dielectric += cylinder(spindleRadius, 2.0 * (size + gap), xyz - vec3(0, 0, dielectricHeight - 4.0));
  dielectric = clamp(dielectric, 0.0, 1.0);
  // This next step is important... we don't want any dielectric material wherever
  // the metal is located, so we just subtract the metal out of the dielectric.
  dielectric -= clamp(metal1 + metal2, 0.0, 1.0);
  
  return vec3(metal1, metal2, dielectric);
}

const float innerRadius = 3.0;
const float wireSize = 0.85; // mm
const float wireGap = 0.15; // mm
const int turnsPerCoil = 122;
const int numCoilPairs = 10;

void mainModel4(out vec4 materials, in vec3 xyz) {
  const float wireSpacing = wireSize + wireGap;
  xyz.z += float(turnsPerCoil) * wireSpacing / 2.0;
  materials.xyz = bifilarElectromagnet(numCoilPairs, innerRadius, wireSize, wireGap, turnsPerCoil, xyz);
}