Skip to content

Instantly share code, notes, and snippets.

@dharmatech
Last active November 10, 2023 01:24
Show Gist options
  • Star 35 You must be signed in to star a gist
  • Fork 2 You must be signed in to fork a gist
  • Save dharmatech/a5e74ef03d98b3ff1c45 to your computer and use it in GitHub Desktop.
Save dharmatech/a5e74ef03d98b3ff1c45 to your computer and use it in GitHub Desktop.

Let's solve the following physics problem using Symbolism, a computer algebra library for C#.

One strategy in a snowball fight is to throw a first snowball at a high angle over level ground. While your opponent is watching the first one, you throw a second one at a low angle and timed to arrive at your opponent before or at the same time as the first one.

Assume both snowballs are thrown with a speed of 25.0 m/s.

The first one is thrown at an angle of 70.0° with respect to the horizontal.

(a) At what angle should the second (lowangle) snowball be thrown if it is to land at the same point as the first?

(b) How many seconds later should the second snowball be thrown if it is to land at the same time as the first?

Define the symbols we'll be using:

var xA = new Symbol("xA");      // position
var yA = new Symbol("yA");

var vxA = new Symbol("vxA");    // velocity
var vyA = new Symbol("vyA");

var vA = new Symbol("vA");

var thA = new Symbol("thA");    // direction

var xB = new Symbol("xB");      
var yB = new Symbol("yB");

var vxB = new Symbol("vxB");    
var vyB = new Symbol("vyB");

var tAB = new Symbol("tAB");    // time traversed
var tAC = new Symbol("tAC");

var ax = new Symbol("ax");      // acceleration
var ay = new Symbol("ay");

var Pi = new Symbol("Pi");

Here, point A is the initial position and B is the peak of trajectory.

The equations we'll be using, based on the standard kinematics equations:

var eqs = new And(

    vxA == vA * cos(thA),
    vyA == vA * sin(thA),

    vxB == vxA + ax * tAB,
    vyB == vyA + ay * tAB,

    xB == xA + vxA * tAB + ax * (tAB ^ 2) / 2,
    yB == yA + vyA * tAB + ay * (tAB ^ 2) / 2
);

Let's define equations for the values that are known:

var vals = new List<Equation>() { xA == 0, yA == 0, vA == 25.0, vyB == 0, ax == 0, ay == -9.8, Pi == Math.PI };

Let's also define a list which is a subset of vals where the value is zero:

var zeros = vals.Where(eq => eq.b == 0).ToList();

OK, let's display our list of equations, eqs:

eqs.DispLong();

Let's re-evaluate the equations with the zeros substituted:

eqs
    .SubstituteEqLs(zeros)
    .DispLong();

OK, now let's use EliminateVariables to eliminate the symbol yB from our set of equations:

eqs
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB)
    .DispLong();

How does EliminateVariables work? It finds the first equation which can be solved for the given symbol. It then solves that particular equation for the symbol and uses this result to substitute all instances of that symbol in the rest of the equations. In this case, there was only one equation with yB.

Now let's eliminate vxA:

eqs
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB, vxA)
    .DispLong();

We can keep eliminating symbols until we've whittled the set of equations down to one:

eqs
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB, vxA, vyA, vxB, tAB)
    .DispLong();

I'd like to solve this equation for thA. However, there are two instances of thA, one in sin and another in cos, and the IsolateVariable method doesn't currently know how to handle that. So let's apply the double-angle formula to combine the sin and cos into a single sin:

eqs
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB, vxA, vyA, vxB, tAB)
    .DeepSelect(DoubleAngleFormulaFunc)
    .DispLong();

Great. Now let's isolate thA:

eqs
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB, vxA, vyA, vxB, tAB)
    .DeepSelect(DoubleAngleFormulaFunc)
    .IsolateVariable(thA)
    .DispLong();

So we see that there are two possible values for thA. Let's save this result and extract the two values of thA:

var expr = eqs
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB, vxA, vyA, vxB, tAB)
    .DeepSelect(DoubleAngleFormulaFunc)
    .IsolateVariable(thA);

var th1 = ((expr as Or).args[0] as Equation).b;
var th2 = ((expr as Or).args[1] as Equation).b;

var th_delta = new Symbol("th_delta");

OK, now we'll add a new equation to eqs. This equation relates th_delta to the difference between the two possible angles. We'll eliminate symbols from this set of equations till we get it down to just one:

eqs
    .Add(th_delta == (th1 - th2).AlgebraicExpand())
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB, vxA, vyA, vxB, tAB)
    .DeepSelect(DoubleAngleFormulaFunc)
    .EliminateVariable(xB)
    .DispLong()

Let's calculate the numeric value:

eqs
    .Add(th_delta == (th1 - th2).AlgebraicExpand())
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB, vxA, vyA, vxB, tAB)
    .DeepSelect(DoubleAngleFormulaFunc)
    .EliminateVariable(xB)                            
    .SubstituteEq(thA == (70).ToRadians())
    .SubstituteEq(Pi == Math.PI)
    .DispLong()

This is approximately -50 degrees. So we know the second snowball should be thrown at 20 degrees.

Now on to part b. Let's eliminate symbols and isolate tAB, the time from point A to point B. We'll save the result in a variable.

var tAB_eq = eqs
    .SubstituteEqLs(zeros)
    .EliminateVariables(yB, vxA, vyA, vxB, xB)
    .IsolateVariable(tAB)
    .DispLong();

We know that thA can be either 70 degrees or 20 degrees. This can be expressed with an Or. We also know that point B is the half-way point so the time to point C is twice the time to B. Combining these with the equation we just found for tAB, we get:

new And(
    new Or(thA == (20).ToRadians(), thA == (70).ToRadians()),
    tAB_eq,
    tAC == tAB * 2)

    .DispLong()

Expanding and eliminating symbols:

new And(
    new Or(thA == (20).ToRadians(), thA == (70).ToRadians()),
    tAB_eq,
    tAC == tAB * 2)
    .LogicalExpand()
    .EliminateVariables(thA, tAB)
    .DispLong()

Let's substitute vals to find the numeric values:

new And(
    new Or(thA == (20).ToRadians(), thA == (70).ToRadians()),
    tAB_eq,
    tAC == tAB * 2)
    .LogicalExpand()
    .EliminateVariables(thA, tAB)
    .SubstituteEqLs(vals)
    .DispLong()

So it takes one of the snowballs about 1.74 seconds to land. It takes the other one 4.79 seconds. If the launcher wants them to get there at the same time, the second one should be thrown 3.05 seconds later.

This example is used as a unit test.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment