forked from gmlewis/irmf
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
Showing
2 changed files
with
368 additions
and
0 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
185 changes: 185 additions & 0 deletions
185
...s/012-bifilar-electromagnet/axial+radial-bifilar-electromagnet-with-solid-core-rot90.irmf
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,185 @@ | ||
/*{ | ||
"author": "Glenn M. Lewis", | ||
"copyright": "Apache-2.0", | ||
"date": "2020-03-05", | ||
"irmf": "1.0", | ||
"materials": ["33CrMoV12","33CrMoV12","dielectric"], | ||
"max": [71,25,25], | ||
"min": [-61,-25,-25], | ||
"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": "axial+radial bifilar electromagnet - full model - metal, dielectric, with solid core, horizontal", | ||
"units": "mm", | ||
"version": "1.0" | ||
}*/ | ||
|
||
#define M_PI 3.1415926535897932384626433832795 | ||
|
||
float coilSquareFace(in float radius, in float size, in float gap, in int nTurns, in float trimStartAngle, in float trimEndAngle, in vec3 xyz) { | ||
// First, trivial reject on the two vertical ends of the coil. | ||
if (xyz.z < -0.5 * size || xyz.z > float(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 >= float(nTurns)) { return 0.0; } | ||
|
||
// For the last two checks, convert angle back to radians. | ||
angle *= (2.0 * M_PI); | ||
|
||
// If the current point is in the first coil, use the trimStartAngle. | ||
if (coilNum < 1.0) { return angle >= trimStartAngle ? 1.0 : 0.0; } | ||
|
||
// If the current point is in the last coil, use the trimEndAngle. | ||
if (coilNum >= float(nTurns) - 1.0 && trimEndAngle > 0.0) { return angle <= trimEndAngle ? 1.0 : 0.0; } | ||
|
||
return 1.0; | ||
} | ||
|
||
float box(vec3 start, vec3 end, float size, in vec3 xyz) { | ||
vec3 ll = min(start, end) - vec3(0.5 * size); | ||
vec3 ur = max(start, end) + vec3(0.5 * size); | ||
if (any(lessThan(xyz, ll))|| any(greaterThan(xyz, ur))) { return 0.0; } | ||
return 1.0; | ||
} | ||
|
||
mat3 rotAxis(vec3 axis, float a) { | ||
// This is from: http://www.neilmendoza.com/glsl-rotation-about-an-arbitrary-axis/ | ||
float s = sin(a); | ||
float c = cos(a); | ||
float oc = 1.0 - c; | ||
vec3 as = axis * s; | ||
mat3 p = mat3(axis.x * axis, axis.y * axis, axis.z * axis); | ||
mat3 q = mat3(c, - as.z, as.y, as.z, c, - as.x, - as.y, as.x, c); | ||
return p * oc + q; | ||
} | ||
|
||
mat4 rotZ(float angle) { | ||
return mat4(rotAxis(vec3(0, 0, 1), angle)); | ||
} | ||
|
||
float wire(vec3 start, vec3 end, float size, in vec3 xyz) { | ||
vec3 v = end - start; | ||
float angle = dot(v, vec3(1, 0, 0)); | ||
xyz -= start; | ||
xyz = (vec4(xyz, 1) * rotZ(angle)).xyz; | ||
return box(vec3(0), vec3(length(v), 0, 0), size, xyz); | ||
} | ||
|
||
float coilPlusConnectorWires(in int wireNum, in int coilNum, in int numCoils, in float inc, in float innerRadius, in float connectorRadius, in float size, in float singleGap, in int nTurns, in vec3 xyz) { | ||
float gap = size + 2.0 * singleGap; | ||
float radiusOffset = float(coilNum - 1); | ||
// next line is experimental - works for numPairs=20: | ||
float spacingAngle = float(numCoils - 4) * inc * atan(2.0 * radiusOffset / float(numCoils - 1)); | ||
mat4 xfm = mat4(1) * rotZ(spacingAngle); | ||
float coilRadius = radiusOffset + innerRadius; | ||
float trimStartAngle = 0.05; // experimental - works somewhat for numPairs=20, but not ideal. | ||
|
||
xyz = (vec4(xyz, 1.0) * xfm).xyz; | ||
float coil = coilSquareFace(coilRadius, size, gap, nTurns, trimStartAngle, 0.0, xyz); | ||
|
||
float tz = float(nTurns) * (size + gap); | ||
float zexit = (float(nTurns) + 10.0) * (size + gap); | ||
|
||
if (coilNum == numCoils) { | ||
// Special case to access the exit wires. | ||
if (wireNum == 1) { | ||
coil += box(vec3(coilRadius, 0, tz), vec3(coilRadius, 0, zexit), size, xyz); | ||
} else { | ||
coil += box(vec3(coilRadius, 0, tz), vec3(connectorRadius, 0, tz), size, xyz); | ||
coil += box(vec3(connectorRadius, 0, - size), vec3(connectorRadius, 0, tz), size, xyz); | ||
coil += box(vec3(connectorRadius, 0, - size), vec3(innerRadius, 0, - size), size, xyz); | ||
coil += box(vec3(innerRadius, 0, - size), vec3(innerRadius, 0, 0), size, xyz); | ||
// Now the other exit wire | ||
coil += box(vec3(-innerRadius, 0, - size), vec3(-innerRadius, 0, 0), size, xyz); | ||
coil += box(vec3(-connectorRadius, 0, - size), vec3(-innerRadius, 0, - size), size, xyz); | ||
coil += box(vec3(-connectorRadius, 0, - size), vec3(-connectorRadius, 0, zexit), size, xyz); | ||
} | ||
} else { | ||
coil += box(vec3(coilRadius, 0, tz), vec3(coilRadius, 0, tz + size), size, xyz); | ||
coil += box(vec3(coilRadius, 0, tz + size), vec3(connectorRadius, 0, tz + size), size, xyz); | ||
coil += box(vec3(connectorRadius, 0, - size), vec3(connectorRadius, 0, tz + size), size, xyz); | ||
coilRadius += size + singleGap; | ||
coil += box(vec3(coilRadius, 0, - size), vec3(connectorRadius, 0, - size), size, xyz); | ||
coil += box(vec3(coilRadius, 0, - size), vec3(coilRadius, 0, 0), size, xyz); | ||
} | ||
|
||
return clamp(coil, 0.0, 1.0); | ||
} | ||
|
||
float cylinder(float radius, float height, in vec3 xyz) { | ||
// First, trivial reject on the two ends of the cylinder. | ||
if (xyz.z < 0.0 || xyz.z > height) { return 0.0; } | ||
|
||
// Then, constrain radius of the cylinder: | ||
float rxy = length(xyz.xy); | ||
if (rxy > radius) { return 0.0; } | ||
|
||
return 1.0; | ||
} | ||
|
||
vec3 arBifilarElectromagnet(int numPairs, float innerRadius, float size, float gap, int numTurns, in vec3 xyz) { | ||
float nTurns = float(numTurns); | ||
float inc = M_PI / float(numPairs); | ||
float connectorRadius = innerRadius + float(numPairs) * (size + gap); | ||
|
||
float metal1 = 0.0; | ||
float metal2 = 0.0; | ||
|
||
mat4 xfm = mat4(1) * rotZ(M_PI); | ||
vec3 xyz180Z = (vec4(xyz, 1.0) * xfm).xyz; | ||
|
||
for(int i = 1; i <= numPairs; i ++ ) { | ||
metal1 += coilPlusConnectorWires(1, i, numPairs, inc, innerRadius, connectorRadius, size, gap, numTurns, xyz); | ||
metal2 += coilPlusConnectorWires(2, i, numPairs, inc, innerRadius, connectorRadius, size, gap, numTurns, xyz180Z); | ||
} | ||
metal1 = clamp(metal1, 0.0, 1.0); | ||
metal2 = clamp(metal2 - metal1, 0.0, 1.0); // Don't double the metals on overlap. | ||
|
||
float dielectricRadius = float(numPairs) * (size + gap) + innerRadius + gap; | ||
float dielectricHeight = 2.0 * (nTurns + 2.0) * (size + gap); | ||
float dielectric = cylinder(dielectricRadius, dielectricHeight, xyz + vec3(0, 0, 2)); | ||
|
||
float solidCore = 2.0 * cylinder(innerRadius - 0.5 * size - gap - 0.05, dielectricHeight, xyz + vec3(0, 0, 2)); | ||
|
||
float spindleRadius = 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 -= 2.0 * cylinder(innerRadius - 0.5 * size - gap, dielectricHeight, xyz + vec3(0, 0, 2)); | ||
// 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(dielectric - metal1 - metal2, 0.0, 1.0); | ||
|
||
return vec3(metal1 + metal2, solidCore, dielectric); | ||
} | ||
|
||
void mainModel4(out vec4 materials, in vec3 xyz) { | ||
xyz.x += 60.0; | ||
materials.xyz = arBifilarElectromagnet(20, 3.0, 0.85, 0.15, 61, xyz.zyx); | ||
} |