Skip to content

Embed URL

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
namespace PSE_eq4_8_eq4_9
{
public static class MathExtensionMethods
{
public static double Sqrt(this double n)
{ return Math.Sqrt(n); }
public static double Sin(this double n)
{ return Math.Sin(n); }
public static double Cos(this double n)
{ return Math.Cos(n); }
public static double Asin(this double n)
{ return Math.Asin(n); }
public static double Acos(this double n)
{ return Math.Acos(n); }
public static double Pow(this double a, double b)
{ return Math.Pow(a, b); }
public static double Degrees(this double n)
{ return 180 * n / Math.PI; }
public static double Radians(this double n)
{ return n * Math.PI / 180; }
}
class Point
{
public double? x, y;
public double X { get { return x.Value; } }
public double Y { get { return y.Value; } }
public override string ToString()
{ return string.Format("Point({0,7:F3}, {1,7:F3})", x, y); }
public static Point FromAngle(double angle, double mag)
{ return new Point(angle.Cos() * mag, angle.Sin() * mag); }
public Point(double? x_val, double? y_val) { x = x_val; y = y_val; }
public static Point operator +(Point a, Point b)
{ return new Point(a.x + b.x, a.y + b.y); }
public static Point operator -(Point a, Point b)
{ return new Point(a.x - b.x, a.y - b.y); }
public static Point operator *(Point p, double n)
{return new Point(p.x * n, p.y * n); }
public static Point operator *(double n, Point p)
{ return new Point(p.x * n, p.y * n); }
public static Point operator /(Point p, double n)
{return new Point(p.x / n, p.y / n); }
public static Point operator /(double n, Point p)
{ return new Point(p.x / n, p.y / n); }
public double Norm()
{ return Math.Sqrt((double)(x * x + y * y)); }
public double ToAngle()
{
return Math.Atan2((double) y, (double) x);
}
}
// two-dimensional motion with constant acceleration
// vf_ = vi_ + a_ t (1)
// (1) vi vi_ = vf_ - a_ t (1.1)
// (1) a a_ = (vf_ - vi_) / t (1.2)
// (1) t t = (vf_ - vi_) / a_ (1.3)
// (1) x component: vfx = vix + ax t (1.x)
// (1) y component: vfy = viy + ay t (1.y)
// (1.x) /. (3) vfx = vi cos(th) + ax t (1.x.1)
// (1.y) /. (4) vfy = vi sin(th) + ay t (1.y.1)
//////////////////////////////////////////////////////////////////////
// rf_ = ri_ + vi_ t + 1/2 a_ t^2 (2)
// (2) vi_ vi_ = (rf_ - ri_ - 1/2 a_ t^2) / t (2.1)
// (2) a_ a_ = (rf_ - ri_ - vi_ t) / (1/2 t^2) (2.2)
//////////////////////////////////////////////////////////////////////
// (2) x components: xf = xi + vxi t + 1/2 ax t^2 (2.3)
//////////////////////////////////////////////////////////////////////
// (2) y components: yf = yi + vyi t + 1/2 ay t^2 (2.4)
// (2.4) t: 0 = 1/2 ay t^2 + vyi t + yi - yf
// t = (- vyi +- Sqrt(vyi^2 - 4 * 1/2 ay (yi - yf))) / (2 * 1/2 ay)
// t = (- vyi +- Sqrt(vyi^2 - 2 ay (yi - yf))) / ay (2.5)
//////////////////////////////////////////////////////////////////////
#region // (2.3) and (2.4) in terms of speed and direction
// vxi = vi cos(th) (3)
// vyi = vi sin(th) (4)
// (2.3) /. (3) xf = xi + vi cos(th) t + 1/2 ax t^2 (2.3.1)
// (2.4) /. (4) yf = yi + vi sin(th) t + 1/2 ay t^2 (2.4.1)
#endregion
//////////////////////////////////////////////////////////////////////
#region // projectile motion - level plane - initial angle or speed known
// projectile motion (ax = 0)
// (2.3.1): xf = xi + vi cos(th) t (2.3.2)
//////////////////////////////////////////////////////////////////////
// projectile motion : level plane (yi = yf)
// (2.4.1): 0 = vi sin(th) t + 1/2 ay t^2
// - vi sin(th) t = 1/2 ay t^2
// - vi sin(th) = 1/2 ay t (2.4.2)
//////////////////////////////////////////////////////////////////////
// (2.3.2) t: t = (xf - xi) / (vi cos(th)) (2.3.2.1)
// (2.4.2) /. (2.3.2.1)
// - vi sin(th) = 1/2 ay (xf - xi) / (vi cos(th))
// - vi sin(th) cos(th) = 1/2 ay (xf - xi) / vi
// - vi^2 sin(th) cos(th) = 1/2 ay (xf - xi)
// vi^2 sin(th) cos(th) = - 1/2 ay (xf - xi)
// vi^2 2 sin(th) cos(th) = - ay (xf - xi)
// vi^2 sin(2 th) = - ay (xf - xi) (2.4.2.1)
//////////////////////////////////////////////////////////////////////
// (2.4.2.1) th: sin(2 th) = - ay (xf - xi) / vi^2
// 2 th = arcsin(- ay (xf - xi) / vi^2)
// th = 1/2 arcsin(- ay (xf - xi) / vi^2) (2.4.2.1.1)
//////////////////////////////////////////////////////////////////////
// (2.4.2.1) vi: vi^2 = - ay (xf - xi) / sin(2 th)
// vi = +- sqrt(- ay (xf - xi) / sin(2 th)) (2.4.2.1.2)
#endregion
//////////////////////////////////////////////////////////////////////
#region // projectile motion incline
// xB = xA + vAx t + 1/2 ax t^2
// yB = yA + vAy t + 1/2 ay t^2
// xB - xA = d cos(phi)
// yB - yA = d sin(phi)
// ax = 0
// xA = 0
// xB = xA + vAx t
// xB - xA = vAx t
// d cos(phi) = vAx t
// cos(phi) = vAx t / d
// cos(phi) / vAx = t / d
// t = d cos(phi) / vAx
// yB - yA = vAy t + 1/2 ay t^2
// d sin(phi) = vAy t + 1/2 ay t^2
// sin(phi) = vAy t / d + 1/2 ay t^2 / d
// sin(phi) = vAy cos(phi) / vAx + 1/2 ay t^2 / d
// sin(phi) = vAy cos(phi) / vAx + 1/2 ay (d cos(phi) / vAx)^2 / d
// sin(phi) = vAy cos(phi) / vAx + 1/2 ay d cos^2(phi) / vAx^2
// d = {sin(phi) - vAy / vAx cos(phi)} / {1/2 ay cos^2(phi) / vAx^2}
#endregion
//////////////////////////////////////////////////////////////////////
class Obj
{
public static bool NotNull(params object[] objects)
{
for (int i = 0; i < objects.Length; i++)
if (objects[i] == null) return false;
return true;
}
static double Discriminant(double a, double b, double c)
{
return b*b - 4*a*c;
}
public static double QuadraticA(double a, double b, double c)
{
return (-b + Math.Sqrt(b * b - 4 * a * c)) / (2 * a);
}
public static double QuadraticB(double a, double b, double c)
{
return (-b - Math.Sqrt(b * b - 4 * a * c)) / (2 * a);
}
public Point position, velocity, acceleration;
public double? time;
public double Time { get { return time.Value; } }
public double? angle; // constraint
public double? speed; // constraint
public double Speed { get { return speed.Value; } }
public double Angle { get { return angle.Value; } }
public override string ToString()
{
return string.Format("Obj(pos={0}, vel={1}, acc={2}, time={3:F3})",
position, velocity, acceleration, time);
}
public Obj AtTime(double t)
{
var obj = new Obj();
obj.time = t;
var dt = (double) (t - time);
obj.acceleration = acceleration;
obj.velocity = velocity + acceleration * dt;
obj.position = position + velocity * dt + 0.5 * acceleration * dt * dt;
return obj;
}
public static Point CalcAcceleration(Obj a, Obj b)
{
var dt = (double)(b.time - a.time);
if (b.velocity != null && a.velocity != null)
return (b.velocity - a.velocity) / dt;
if (b.position != null && a.position != null && a.velocity != null)
return (b.position - a.position - a.velocity * dt) / (1 / 2 * dt * dt);
throw new Exception();
}
public static double CalcTime(Obj a, Obj b, int solution=0)
{
if (NotNull(a.position.x, b.position.x, a.velocity.x) && a.acceleration.x == 0.0)
return (b.position.X - a.position.X) / a.velocity.X;
if (b.velocity.x != null && a.velocity.x != null && a.acceleration.x != null && a.acceleration.x != 0.0)
return (double) ((b.velocity.x - a.velocity.x) / a.acceleration.x);
if (b.velocity.y != null && a.velocity.y != null && a.acceleration.y != null && a.acceleration.y != 0.0)
return (double) ((b.velocity.y - a.velocity.y) / a.acceleration.y);
Func<double, double, double, double> quadratic;
if (solution == 0)
quadratic = QuadraticA;
else
quadratic = QuadraticB;
if (a.position.x != null &&
b.position.x != null &&
a.velocity.x != null &&
a.acceleration.x != null &&
a.acceleration.x != 0.0)
return
quadratic(
(double)(0.5 * a.acceleration.x),
(double)a.velocity.x,
(double)(a.position.x - b.position.x));
if (a.position.y != null &&
b.position.y != null &&
a.velocity.y != null &&
a.acceleration.y != null &&
a.acceleration.y != 0.0)
return
quadratic(
(double) (0.5 * a.acceleration.y),
(double) a.velocity.y,
(double) (a.position.y - b.position.y));
throw new Exception();
}
public static List<double> CalcTimes(Obj a, Obj b)
{
var result = new List<double>();
if (a.position.y != null &&
b.position.y != null &&
a.velocity.y != null &&
a.acceleration.y != null)
{
var discriminant =
Discriminant(
(double)(0.5 * a.acceleration.y),
(double)a.velocity.y,
(double)(a.position.y - b.position.y));
if (discriminant > 0)
{
result.Add(
QuadraticA(
(double)(0.5 * a.acceleration.y),
(double)a.velocity.y,
(double)(a.position.y - b.position.y)));
result.Add(
QuadraticB(
(double)(0.5 * a.acceleration.y),
(double)a.velocity.y,
(double)(a.position.y - b.position.y)));
return result;
}
if (discriminant == 0.0)
{
result.Add(
QuadraticA(
(double)(0.5 * a.acceleration.y),
(double)a.velocity.y,
(double)(a.position.y - b.position.y)));
return result;
}
throw new Exception("no real values");
}
throw new Exception();
}
public static Point CalcInitialVelocity(Obj a, Obj b, int solution=0)
{
if (NotNull(a.time, b.time))
{
var dt = (double)(b.time - a.time);
if (NotNull(b.velocity.x, b.velocity.y, b.acceleration.x, b.acceleration.y))
return b.velocity - b.acceleration * dt;
if (NotNull(b.position.x, b.position.y, b.acceleration.x, b.acceleration.y))
return (b.position - a.position - 0.5 * b.acceleration * dt * dt) / dt;
}
if (a.acceleration.x == 0 &&
b.acceleration.x == 0 &&
a.position.y == b.position.y)
{
if (a.speed != null)
{
a.angle = 0.5 * Math.Asin(-a.acceleration.Y * (b.position.X - a.position.X) / a.Speed.Pow(2));
return Point.FromAngle(a.Angle, a.Speed);
}
if (a.angle != null)
{
a.speed = Math.Sqrt(-a.acceleration.Y * (b.position.X - a.position.X) / Math.Sin(2 * a.Angle));
return Point.FromAngle(a.Angle, a.Speed);
}
}
if (a.position.x == 0 &&
a.position.y == 0 &&
b.position.x.HasValue &&
b.position.y.HasValue &&
a.acceleration.x == 0 &&
b.acceleration.x == 0 &&
a.acceleration.y.HasValue)
{
if (a.speed.HasValue)
{
var xf = b.position.X;
var yf = b.position.Y;
var vi = a.Speed;
var ay = a.acceleration.Y;
if (solution == 0)
a.angle = (Math.Asin((yf - ay * xf.Pow(2) / vi.Pow(2)) / Math.Sqrt(xf.Pow(2) + yf.Pow(2))) + Math.Atan(yf / xf)) / 2.0;
if (solution == 1)
a.angle = (Math.PI - Math.Asin((yf - ay * xf.Pow(2) / vi.Pow(2)) / Math.Sqrt(xf.Pow(2) + yf.Pow(2))) + Math.Atan(yf / xf)) / 2.0;
return Point.FromAngle(a.Angle, a.Speed);
}
}
throw new Exception();
}
//public double CalcRange()
//{
// return
// Math.Pow(velocity.Norm(), 2) *
// Math.Sin(2 * velocity.ToAngle()) /
// Math.Abs((double)acceleration.y);
//}
public Point ProjectileInclineIntersection(double phi)
{
if (NotNull(phi, velocity.x, velocity.y, acceleration.y)
&&
acceleration.x == 0.0)
{
var d =
(phi.Sin() - velocity.y / velocity.x * phi.Cos())
/
(0.5 * acceleration.y * phi.Cos().Pow(2) / ((double)velocity.x).Pow(2));
return
new Point(
position.x + d * phi.Cos(),
position.y + d * phi.Sin());
}
throw new Exception();
}
// ar = v^2 / r
public double CalcRadialAcceleration(double radius)
{ return velocity.Norm().Pow(2) / radius; }
}
class Program
{
// A placekicker must kick a football from a point 36.0 m
// (about 40 yards) from the goal, and half the crowd
// hopes the ball will clear the crossbar, which is 3.05 m
// high. When kicked, the ball leaves the ground with a
// speed of 20.0 m/s at an angle of 53.0° to the horizontal.
//
// (a) By how much does the ball clear or fall short of
// clearing the crossbar?
//
// (b) Does the ball approach the crossbar while
// still rising or while falling?
// PSE PRO 4.19
// xi = 0
// yi = 0
// vix = 20 cos(53 deg)
// viy = 20 sin(53 deg)
// ax = 0
// ay = -9.8
// xf = 36
// yf = ?
// t = ?
// (2.3) xf = xi + vxi t + 1/2 ax t^2
// ax=0 xf = xi + vxi t
// t = (xf - xi) / vxi
static void Main(string[] args)
{
var obj_A = new Obj();
obj_A.position = new Point(0, 0);
obj_A.velocity = Point.FromAngle((53.0).Radians(), 20);
obj_A.acceleration = new Point(0, -9.8);
obj_A.time = 0;
var obj_B = new Obj();
obj_B.position = new Point(36, null);
obj_B.velocity = new Point(obj_A.velocity.x, null);
obj_B.acceleration = new Point(0, -9.8);
obj_B.time = null;
obj_B.time = Obj.CalcTime(obj_A, obj_B);
obj_B = obj_A.AtTime(obj_B.Time);
Console.WriteLine("object at A: " + obj_A);
Console.WriteLine("object at B: " + obj_B);
var field_goal_height = 3.05;
Console.WriteLine("ball clears goal by {0:F3} meters",
(obj_B.position.y - field_goal_height));
Console.WriteLine("ball approaches the goal while {0}",
obj_B.velocity.y > 0 ? "rising" : "falling");
Console.ReadLine();
}
}
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Something went wrong with that request. Please try again.