// Predefined design scripts, selectable from the interface. Each is a
// MATLAB-syntax script that fills in the `design` struct; see design.ts for
// the fields. A script is self-contained: the shape functions it uses are
// ordinary MATLAB functions written at the bottom of the script itself.

export interface Preset {
  id: string;
  label: string;
  source: string;
}

const flatRectangle = `% A flat rectangular sheet of stainless steel, hanging freely.
% The script is an ordinary MATLAB program that fills in a struct named
% design. Lengths are meters; the constants at the top are the only units.
% A field can be a number or a function handle @(x, y, r, theta) of the
% point, measured from the center of the bounding box; the handles are
% applied to whole grid matrices, so use elementwise operators (.*, ./, .^)
% inside them. Anything else a script needs it defines itself, as MATLAB
% functions at the bottom; the other presets show how. The decay of the
% sound is not prescribed: it follows from the loss factor of the metal, its
% thermal constants (thermoelastic damping), the radiation of each mode into
% the air.
inch = 0.0254;                  % m
mm = 1e-3;                      % m

design.name = 'Flat steel rectangle';
design.material.E = 170e9;      % Young's modulus [Pa]
design.material.nu = 0.29;      % Poisson ratio
design.material.rho = 7900;     % density [kg/m^3]
design.material.lossFactor = 3e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 17e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 500;            % specific heat [J/(kg K)]  }
design.material.k = 16;              % conductivity [W/(m K)]    }
design.bounds.width = 12*inch;  % bounding box [m]
design.bounds.height = 18*inch;
design.thickness = 0.79*mm;
`;

const aluminumSquare = `% A flat 12 inch square of 6061-T6 aluminum, 1.02 mm thick.
% Aluminum is a third the density of steel and a third as stiff, so
% the bending frequencies land near those of a steel plate of the same size
% and thickness. Its high thermal expansion and conductivity give it the
% strongest thermoelastic damping of the three metals here.
inch = 0.0254; mm = 1e-3;
W = 12*inch;

design.name = 'Aluminum 6061-T6 square';
design.material.E = 68.9e9;     % Young's modulus [Pa]
design.material.nu = 0.33;      % Poisson ratio
design.material.rho = 2700;     % density [kg/m^3]
design.material.lossFactor = 1e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 23.6e-6;     % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 896;            % specific heat [J/(kg K)]  }
design.material.k = 167;             % conductivity [W/(m K)]    }
design.bounds.width = W;
design.bounds.height = W;
design.thickness = 1.02*mm;
`;

const bowedRectangle = `% The same steel sheet with a 3 mm parabolic bow across its width, the
% kind of bow that rolling and handling leave in real sheet metal. The bow
% makes the sheet a shallow shell: the linear modes couple bending to
% stretching, and the full simulation oil-cans and wobbles when struck.
inch = 0.0254; mm = 1e-3;
W = 12*inch; H = 18*inch;

design.name = 'Bowed steel rectangle';
design.material.E = 170e9;
design.material.nu = 0.29;
design.material.rho = 7900;
design.material.lossFactor = 3e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 17e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 500;            % specific heat [J/(kg K)]  }
design.material.k = 16;              % conductivity [W/(m K)]    }
design.bounds.width = W;
design.bounds.height = H;
design.thickness = 0.79*mm;
design.profile = @(x, y, r, theta) bow(x, W, 3*mm, 0.3);   % stress-free height [m]

function v = bow(s, span, sag, offset)
% A parabola of the given sag across a span centered on 0; s runs along the
% span, and the offset in -1..1 moves the crest toward one end.
t = 2*s/span;
v = sag * (1 - t.^2) .* (1 + offset*t);
end
`;

const flatDisc = `% A flat 16 inch bronze disc, hanging freely.
inch = 0.0254; mm = 1e-3;
D = 16*inch;

design.name = 'Flat bronze disc';
design.material.E = 110e9;      % bronze
design.material.nu = 0.34;
design.material.rho = 8800;
design.material.lossFactor = 2e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 18e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 370;            % specific heat [J/(kg K)]  }
design.material.k = 40;              % conductivity [W/(m K)]    }
design.bounds.width = D;
design.bounds.height = D;
design.inside = @(x, y, r, theta) r <= D/2;
design.thickness = 1.0*mm;
`;

const crash16 = `% A 16 inch flat crash: bronze, lathed to a taper from 1.2 mm at the
% center to 0.7 mm at the edge. Flat, so the linear modes are pure bending
% modes; the tapered profile is what makes it a crash rather than a ride.
inch = 0.0254; mm = 1e-3;
D = 16*inch; R = D/2;
tCenter = 1.2*mm; tEdge = 0.7*mm;

design.name = '16 inch flat crash';
design.material.E = 110e9;
design.material.nu = 0.34;
design.material.rho = 8800;
design.material.lossFactor = 2e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 18e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 370;            % specific heat [J/(kg K)]  }
design.material.k = 40;              % conductivity [W/(m K)]    }
design.bounds.width = D;
design.bounds.height = D;
design.inside = @(x, y, r, theta) r <= R;
design.thickness = @(x, y, r, theta) tCenter + (tEdge - tCenter) * ramp(40*mm, R, r);

function v = ramp(a, b, x)
% 0 at x = a, 1 at x = b, clamped to that range outside.
v = min(max((x - a) ./ (b - a), 0), 1);
end
`;

const ride20 = `% A 20 inch flat ride: thicker bronze, 2.2 mm at the center to 1.4 mm at
% the edge. The thickness raises every frequency in proportion, and the
% ride has a clearer ping and a longer decay than the crash.
inch = 0.0254; mm = 1e-3;
D = 20*inch; R = D/2;

design.name = '20 inch flat ride';
design.material.E = 110e9;
design.material.nu = 0.34;
design.material.rho = 8800;
design.material.lossFactor = 2e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 18e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 370;            % specific heat [J/(kg K)]  }
design.material.k = 40;              % conductivity [W/(m K)]    }
design.bounds.width = D;
design.bounds.height = D;
design.inside = @(x, y, r, theta) r <= R;
design.thickness = @(x, y, r, theta) 2.2*mm + (1.4*mm - 2.2*mm) * ramp(50*mm, R, r);

function v = ramp(a, b, x)
v = min(max((x - a) ./ (b - a), 0), 1);
end
`;

const crashBell16 = `% A 16 inch crash with a shallow bow and a cosine bell. The
% profile makes the cymbal a shallow shell: the linear model couples the
% bending to the stretching (a larger eigenproblem, so the modes take
% longer), and the full simulation shows the asymmetric nonlinearity of a
% curved plate. The slopes are kept small; the interface reports the
% largest one.
inch = 0.0254; mm = 1e-3;
D = 16*inch; R = D/2;
bellR = 55*mm; bellH = 8*mm; bowH = 10*mm;
tCenter = 1.2*mm; tEdge = 0.7*mm;

design.name = '16 inch crash with bell';
design.material.E = 110e9;
design.material.nu = 0.34;
design.material.rho = 8800;
design.material.lossFactor = 2e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 18e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 370;            % specific heat [J/(kg K)]  }
design.material.k = 40;              % conductivity [W/(m K)]    }
design.bounds.width = D;
design.bounds.height = D;
design.inside = @(x, y, r, theta) r <= R;
design.thickness = @(x, y, r, theta) tCenter + (tEdge - tCenter) * ramp(bellR, R, r);
design.profile = @(x, y, r, theta) bowH * (1 - (r/R).^2) + dome(r, bellR, bellH);

function v = ramp(a, b, x)
v = min(max((x - a) ./ (b - a), 0), 1);
end

function v = dome(d, radius, height)
% A cosine dome of the given height and base radius, with zero slope both at
% its rim and at its apex; d is the distance from its center.
v = 0.5 * height * (1 + cos(pi * d / radius)) .* (d < radius);
end
`;

const crashHammered16 = `% The 16 inch crash with bell, hammered: 160 shallow dents of 0.4 mm on
% the bow, laid out on a sunflower spiral so that no two align, or, with
% randomDents set to 1, scattered at random over the bow (dentSeed picks the
% pattern, and the same seed gives the same cymbal). A hammered cymbal is
% no longer axisymmetric, so its degenerate mode pairs split and its modes
% mix, which is part of what makes a cymbal sound like one. The dents here
% are what a script can express, a change of the rest shape; the residual
% stresses of real hammering are not modeled.
inch = 0.0254; mm = 1e-3;
D = 16*inch; R = D/2;
bellR = 55*mm; bellH = 8*mm; bowH = 10*mm;
tCenter = 1.2*mm; tEdge = 0.7*mm;
nDents = 160; dentR = 12*mm; dentDepth = 0.4*mm;
randomDents = 1;                % 0: the spiral; 1: random places, overlapping as hammer marks do
dentSeed = 1;                   % any positive integer; changes the random pattern

design.name = '16 inch hammered crash';
design.material.E = 110e9;
design.material.nu = 0.34;
design.material.rho = 8800;
design.material.lossFactor = 2e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 18e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 370;            % specific heat [J/(kg K)]  }
design.material.k = 40;              % conductivity [W/(m K)]    }
design.bounds.width = D;
design.bounds.height = D;
design.inside = @(x, y, r, theta) r <= R;
design.thickness = @(x, y, r, theta) tCenter + (tEdge - tCenter) * ramp(bellR, R, r);
design.profile = @(x, y, r, theta) bowH * (1 - (r/R).^2) + dome(r, bellR, bellH) + dents(x, y, bellR + dentR, R - dentR, nDents, dentR, dentDepth, randomDents, dentSeed);

function v = ramp(a, b, x)
v = min(max((x - a) ./ (b - a), 0), 1);
end

function v = dome(d, radius, height)
% A cosine dome of the given height and base radius, with zero slope both at
% its rim and at its apex; d is the distance from its center.
v = 0.5 * height * (1 + cos(pi * d / radius)) .* (d < radius);
end

function v = dents(x, y, r0, r1, n, rd, depth, random, seed)
% n hammer marks between the radii r0 and r1: Gaussian dents of 1/e radius
% rd, placed on a sunflower spiral (equal area per mark, golden-angle
% spacing) so that the pattern has no symmetry, or at random with uniform
% density over the annulus, from a small congruential generator so that a
% seed always gives the same pattern.
v = 0 * x;
ga = pi * (3 - sqrt(5));
s = seed;
for k = 1:n
  if random
    s = mod(s * 1664525 + 1013904223, 4294967296); u = s / 4294967296;
    s = mod(s * 1664525 + 1013904223, 4294967296); w = s / 4294967296;
    rk = sqrt(r0^2 + (r1^2 - r0^2) * u);
    tk = 2 * pi * w;
  else
    rk = sqrt(r0^2 + (r1^2 - r0^2) * (k - 0.5) / n);
    tk = k * ga;
  end
  dx = x - rk * cos(tk); dy = y - rk * sin(tk);
  v = v - depth * exp(-(dx.^2 + dy.^2) / rd^2);
end
end
`;

const ellipse = `% An elliptical steel plate. The outline is any test on the point;
% an ellipse inscribed in the bounding box here. The degenerate mode pairs
% of a disc split on an ellipse, which the mode browser shows.
inch = 0.0254; mm = 1e-3;
A = 14*inch; B = 10*inch;

design.name = 'Steel ellipse';
design.material.E = 170e9;
design.material.nu = 0.29;
design.material.rho = 7900;
design.material.lossFactor = 3e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 17e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 500;            % specific heat [J/(kg K)]  }
design.material.k = 16;              % conductivity [W/(m K)]    }
design.bounds.width = A;
design.bounds.height = B;
design.inside = @(x, y, r, theta) (x/(A/2)).^2 + (y/(B/2)).^2 <= 1;
design.thickness = 1.0*mm;
`;

const thunderSheet = `% A thunder sheet: a large thin steel sheet, hung from a corner in reality.
% Its modes are dense and low, and the full simulation turns a hard hit into
% a roar.
% Long at the default spacing; a coarser grid keeps it interactive.
inch = 0.0254; mm = 1e-3;

design.name = 'Thunder sheet';
design.material.E = 200e9;
design.material.nu = 0.3;
design.material.rho = 7850;
design.material.lossFactor = 3e-4;   % internal friction of the metal (typical published value)
design.material.alpha = 12e-6;       % thermal expansion [1/K]   } for the thermoelastic damping
design.material.cp = 460;            % specific heat [J/(kg K)]  }
design.material.k = 50;              % conductivity [W/(m K)]    }
design.bounds.width = 24*inch;
design.bounds.height = 36*inch;
design.thickness = 0.5*mm;
`;

export const PRESETS: Preset[] = [
  { id: "flat-rectangle", label: "Flat steel rectangle", source: flatRectangle },
  { id: "aluminum-square", label: "Aluminum 6061-T6 square", source: aluminumSquare },
  { id: "bowed-rectangle", label: "Bowed steel rectangle", source: bowedRectangle },
  { id: "flat-disc", label: "Flat bronze disc", source: flatDisc },
  { id: "crash-16", label: "16 inch flat crash", source: crash16 },
  { id: "ride-20", label: "20 inch flat ride", source: ride20 },
  { id: "crash-bell-16", label: "16 inch crash with bell", source: crashBell16 },
  { id: "crash-hammered-16", label: "16 inch hammered crash", source: crashHammered16 },
  { id: "ellipse", label: "Steel ellipse", source: ellipse },
  { id: "thunder-sheet", label: "Thunder sheet", source: thunderSheet },
];

export const DEFAULT_PRESET = "crash-16";

export function presetById(id: string): Preset | undefined {
  return PRESETS.find((p) => p.id === id);
}