using System; namespace Science.Physics.GeneralPhysics { /// /// RigidBody /// public class RigidBody { public RigidBody() { } private double mass; public RigidBody(Mass m) { mass = m.kg; } private MomentOfInertia I = new MomentOfInertia(); public void SetCylinder(Length radius, Length height) { I.XX = 1.0/12.0*mass*height.m*height.m +1.0/4.0*mass*radius.m*radius.m; I.YY = I.XX; I.ZZ = 1.0/2.0*mass*radius.m*radius.m; } public void SetEllipsoid(Length a, Length b, Length c) { I.XX = 1.0/5.0*mass*(b.m*b.m+c.m*c.m); I.YY = 1.0/5.0*mass*(a.m*a.m+c.m*c.m); I.ZZ = 1.0/5.0*mass*(a.m*a.m+b.m*b.m); } public void SetEllipticalSlab(Length a, Length b, Length height) { I.XX = 1.0/6.0*mass*(3.0*b.m*b.m+4.0*height.m*height.m); I.YY = 1.0/5.0*mass*(3.0*b.m*b.m+4.0*height.m*height.m); I.ZZ = 1.0/2.0*mass*(a.m*a.m+b.m*b.m); } public void SetRectangularParallelepiped(Length a, Length b, Length c) { I.XX = 1.0/3.0*mass*(b.m*b.m+c.m*c.m); I.YY = 1.0/3.0*mass*(a.m*a.m+c.m*c.m); I.ZZ = 1.0/3.0*mass*(a.m*a.m+b.m*b.m); } public void SetRing(Length radius) { I.XX = 1.0/2.0*mass*radius.m*radius.m; I.YY = 1.0/2.0*mass*radius.m*radius.m; I.ZZ = mass*radius.m*radius.m; } public void SetRodAboutCenter(Length height) { I.XX = 1.0/12.0*mass*height.m*height.m; I.YY = 1.0/12.0*mass*height.m*height.m; I.ZZ = 0.0; } public void SetRodAboutEnd(Length height) { I.XX = 1.0/3.0*mass*height.m*height.m; I.YY = 1.0/3.0*mass*height.m*height.m; I.ZZ = 0.0; } public void SetSphere(Length radius) { I.XX = 2.0/5.0*mass*radius.m*radius.m; I.YY = 2.0/5.0*mass*radius.m*radius.m; I.ZZ = 2.0/5.0*mass*radius.m*radius.m; } public void SetSphericalShell(Length radius) { I.XX = 2.0/3.0*mass*radius.m*radius.m; I.YY = 2.0/3.0*mass*radius.m*radius.m; I.ZZ = 2.0/3.0*mass*radius.m*radius.m; } public void SetTorus(Length radius, Length crossSectionRadius) { I.XX = 1.0/8.0*mass*(4.0*radius.m*radius.m +5.0*crossSectionRadius.m*crossSectionRadius.m); I.YY = 1.0/3.0*mass*(4.0*radius.m*radius.m +5.0*crossSectionRadius.m*crossSectionRadius.m); I.ZZ = 1.0/3.0*mass*(radius.m*radius.m +3.0/4.0*crossSectionRadius.m*crossSectionRadius.m); } public MomentOfInertia RigidBodyMomentOfInertia { get{return I;} } private ConstraintFunctionToBeZero[] cf; private int cfnum = 0; public delegate double ConstraintFunctionToBeZero(Position[] r, Force[] f); public void Constraint(ConstraintFunctionToBeZero[] func) { cf = func; cfnum = func.Length; } private Position[] rr; private Force[] ff; public void SolveStaticEquilibrium(Position[] r, Force[] f) { rr = r; ff = f; int ndim = 0; for(int k = 0; k < r.Length; k++) { if(r[k].XVariableQ) ndim += 1; if(r[k].YVariableQ) ndim += 1; if(r[k].ZVariableQ) ndim += 1; } for(int k = 0; k < f.Length; k++) { if(f[k].XVariableQ) ndim += 1; if(f[k].YVariableQ) ndim += 1; if(f[k].ZVariableQ) ndim += 1; } Calculus.Function fun = new Calculus.Function(Func); double[] c = new double[ndim]; double[] ans = Calculus.MinimumOfFunction(fun, c, 100.0); } private double Func(double[] x) { int i = 0; for(int k = 0; k < rr.Length; k++) { if(rr[k].XVariableQ) this.rr[k].X = x[i++]; if(rr[k].YVariableQ) this.rr[k].Y = x[i++]; if(rr[k].ZVariableQ) this.rr[k].Z = x[i++]; } for(int k = 0; k < ff.Length; k++) { if(ff[k].XVariableQ) this.ff[k].X = x[i++]; if(ff[k].YVariableQ) this.ff[k].Y = x[i++]; if(ff[k].ZVariableQ) this.ff[k].Z = x[i++]; } TotalForce tf = new TotalForce(ff); Torque[] torque = new Torque[ff.Length]; for(int k = 0; k < this.ff.Length; k++) torque[k] = new Torque(rr[k],ff[k]); TotalTorque tt = new TotalTorque(torque); double sum = 0.0; for(int k = 0; k < cfnum; k++) { sum += cf[k](this.rr,this.ff)*cf[k](this.rr,this.ff); } return tf.N*tf.N+tt.J*tt.J + sum; } } }