TODO: Automatically refine the mesh
Closed this issue · 3 comments
bangerth commented
TODO: Automatically refine the mesh
bangerth commented
From Jason:
double c1 = jThinPlateBendingSpeedOfSound(freq, h, rho, E.r(), v);
double c2 = sqrt(T.r() / (rho*h)); // speed of sound in stretched membrane (m/s)
// T=tension (N/m), rho-density(kg/m^3), h=thickness (m)
double c;
if (c1 > c2)
c = c1;
else
c = c2;
double wavelength = c / freq;
#define MIN_DISCRETIZATION_rect_fea_d 15
// minimum # of discretization in x or y direction for this solver
#define DOF_MAX_rect_fea_d 1000
// # maximum of degrees of freedom for this solver
int Nx = (int)(10.0 * Lx / wavelength); // # of discretizations in x direction
int Ny = (int)(10.0 * Ly / wavelength); // # of discretizations in y direction
if (Nx < MIN_DISCRETIZATION_rect_fea_d)
Nx = MIN_DISCRETIZATION_rect_fea_d;
if (Ny < MIN_DISCRETIZATION_rect_fea_d)
Ny = MIN_DISCRETIZATION_rect_fea_d;
if (Nx*Ny>DOF_MAX_rect_fea_d)// ok, we've got to reduce the size
{
double a = Nx*Ny;
double s = sqrt(a / DOF_MAX_rect_fea_d);
Nx = (int)(Nx / s);
Ny = (int)(Ny / s);
if (Nx < MIN_DISCRETIZATION_rect_fea_d)
{
Nx = MIN_DISCRETIZATION_rect_fea_d;
if (Nx*Ny>DOF_MAX_rect_fea_d)
Ny = (int)((double)DOF_MAX_rect_fea_d / (double)Nx);
}
if (Ny < MIN_DISCRETIZATION_rect_fea_d)
{
Ny = MIN_DISCRETIZATION_rect_fea_d;
if (Nx*Ny>DOF_MAX_rect_fea_d)
Nx = (int)((double)DOF_MAX_rect_fea_d / (double)Ny);
}
}
// computes the speed of sound of a shear or vending wave in a thin plate.
// send:
// freq frequency (Hz)
// h heigth or thickness of sheet (m)
// rho density of material (kg/m^3)
// E youngs modulus of material (Pa)
// v Poissons ratio of material (dimensionless)
// (I think v=0 for a sheet of woven material or felt is correct)
// returns:
// speed of sound in m/s
double jThinPlateBendingSpeedOfSound(double freq, double h, double rho, double E, double v)
{
double I_w, w2,a;
I_w = h*h*h / 12; // I/width
// from http://www.engineeringtoolbox.com/area-moment-inertia-d_1328.html
w2 = two_pi*freq;
w2 *= w2;
// from: http://scholar.lib.vt.edu/theses/available/etd-09172001-144121/unrestricted/Blanc_MsThesis.pdf
// speed of sound = (EI*w^2/rho*h*width)^(1/4)
a = E*I_w*w2 / (rho*h);
a = pow(a, 0.25);
return a;
}
bangerth commented
For the record, the formula to compute a
in jThinPlateBendingSpeedOfSound
is missing the factor 1/(1-nu^2)
in the computation of E*I_w
when describing the D
factor in the equations.