CITE THIS NOTEBOOK: Using vector variables and changing variables in ODEs by Peter Burbery. Wolfram Community FEB 4 2023.
Mathematica now directly supports vector variables in ODES. This capability was introduced in version 13.0. I will use this for kinematics.
Exercise 9
Exercise 9
This is from the book A First Course in Differential Equations by J. David Logan, Third Edition published by Springer in 2011 as part of their series Undergraduate texts in Mathematics.
The differential equation can be converted into an antiderivative problem for a new dependent variable using the substitution . Find the differential equation for , solve it, and then determine the general solution of the original equation.
x'=3x+
-t
y=y(t)
x(t)=y(t)
3t
y
x=x(t)
Source: Chapter 1 First-Order Differential Equations Section 1.2 Antiderivatives Exercise 9 page 21
I will solve the problem with DSolve and DSolveValue first.
In[]:=
DSolve[{x'[t]==3x[t]+},x[t],t]
-t
In[]:=
DSolveValue[{x'[t]==3x[t]+},x[t],t]
-t
We determined the general solution of the original equation.
x=x(t)
I figured out how to show the substitution with a new function introduced in 13.1, DSolveChangeVariables. Version 13.1 introduced DSolveChangeVariables for changing variables in differential equations. For more information, Stephen Wolfram describes it in the release announcement for 13.1 section College Calculus
In[]:=
DSolveChangeVariables[Inactive[DSolve][{x'[t]==3x[t]+},x[t],t],y,t,x[t]==y[t]]
-t
3t
In[]:=
DSolveChangeVariables[Inactive[DSolve][{x'[t]==3x[t]+},x[t],t],y(*thenewdependentvariableisy,thepreviousdependentvariablewasx.*),t(*Theindependentvariableistfortime.*),(*Thisisthesubstitution*)x[t]==y[t]]
-t
3t
In[]:=
Activate
We apply the reverse substitution.
In[]:=
y[t]/.Activate
3t
CodeEquivalentQy[t]/.Activate,
3t
The original solution we obtained directly was this.
In[]:=
x[t]/.DSolve[{x'[t]==3x[t]+},x[t],t]
-t
I install the paclet CodeEquivalenceUtilities and load the context with Needs.
In[]:=
PacletInstall["Wolfram/CodeEquivalenceUtilities"]
In[]:=
Needs["Wolfram`CodeEquivalenceUtilities`"]
I will check if the direct solution and the solution with substitution are code equivalent.
In[]:=
CodeEquivalentQx[t]/.DSolve[{x'[t]==3x[t]+},x[t],t](*thedirectsolution*),y[t]/.Activate(*thesubstitutionsolution*)
-t
3t
I like to see the details with EquivalenceTestData.
In[]:=
EquivalenceTestDatax[t]/.DSolve[{x'[t]==3x[t]+},x[t],t](*thedirectsolution*),y[t]/.Activate(*thesubstitutionsolution*)
-t
3t
I will do another exercise.
Exercise 12
Exercise 12
This is exercise 12 from the same page and the same section and the same chapter and the same book as exercise 9.
An enemy cannon at horizontal distance from a fort can fire a cannon ball from the top of a hill at height above the ground level with a muzzle velocity . How high should the wall of the fort be to guarantee that a cannon ball does not go over the wall? Observe that the enemy can adjust the angle of its shot. Hint: Ignoring air resistance, the governing equating follow from resolving Newton's second law for the horizontal and vertical components of the force: and .
L
H
v
θ
mx''=0
my''=-mg
I will use vector variables. Stephen Wolfram writes about the introduction of this in the section As Well As Lots of Other Math… of his release announcement for Version 13.
scratch
scratch
working
working
I will use as a vector with 2 components. The first component will be the horizontal direction as a function of time , . The second component will be the vertical direction as a function of time , . We have is a vector variable.
r(t)
x
t
x(t)
y
t
y(t)
r(t)
r(t){x(t),y(t)}
We need to put the components together to form a vector.
First we have an ODE for the horizontal component of the force.
m(t)0
′′
x
Next, we have an ODE for the vertical component of the force. This time the force downward is given by the magnitude gravitational potential field averaged across the Earth's surface, times the mass, .
g
m
F-mg
We can phrase this as a differential equation.
m(t)g(-m)
′′
y
We have a system of ODEs that we will build a vector with.
m(t)0m(t)g(-m)
′′
x
′′
y
We have
m ′′ x |
m ′′ y |
0 |
-mg |
We can rewrite this by defining r.
r(t)
x(t) |
y(t) |
Now
m(t)
′′
r
m ′′ x |
m ′′ y |
0 |
-mg |
Vector ODE:
m(t)
′′
r
0 |
-mg |
Solve the vector ODE:
In[]:=
DSolveValue[{mr''[t]=={0,-mg}},r[t]∈Vectors[2],t]
I wanted to check this answer is correct so I consulted the textbook Physics: Principles with Applications by Douglas C. Giancoli, 7th edition published by Pearson Education Inc. in 2014.
Source: Chapter 3: Kinematics in Two Dimensions; Vectors Section 3-6 Solving Projectile Motion Problems page 60 Table 3-2 Kinematic Equations for Projectile Motion ( positive upward; =0, =-g=-9.80m/)
y
a
x
a
y
2
s
These formulas are from the book.
The variable represents the initial horizontal position at time . The variable represents the initial horizontal velocity at time . (t).
x
0
t=0
v
x0
t=0
v
x0
′
x
x+t
x
0
v
x0
Source: the column in Table 3-2 mentioned above named Horizontal Motion (=0, =constant) Equation 2-11b
a
x
v
x
The variable represents the initial vertical position at time . The variable represents the initial vertical velocity at time . (t).
y
0
t=0
v
y0
t=0
v
y0
′
y
y+t-g
y
0
v
y0
1
2
2
t
Source: the column in Table 3-2 mentioned above named Vertical (=-g=constant) Equation 2-11b with the note for † If is taken positive downward, the minus (-) signs in front of g become + signs.
†
Motion
a
y
y
I can add the initial conditions for velocity.
′
r
v x |
v y |
′
r
v x0 |
v y0 |
Add the initial condition for position.
r(t)
x(t) |
y(t) |
r(0)
x 0 |
y 0 |
In[]:=
DSolveValue[{mr''[t]=={0,-mg},r'[0]=={,},r[0]=={,}},r[t]∈Vectors[2],t]
v
x0
v
y0
x
0
y
0
Is this the same as the formulas from the book? CodeEquivalentQ knows whether or not it is.
In[]:=
CodeEquivalentQ+t,+t-g(*Giancoli'sequations*),DSolveValue[{mr''[t]=={0,-mg},r'[0]=={,},r[0]=={,}},r[t]∈Vectors[2],t]
x
0
v
x0
y
0
v
y0
1
2
2
t
v
x0
v
y0
x
0
y
0
The derivation is correct.
The equation for parabolic motion is
y==x-
v
y0
v
x0
g
2
2
v
x0
2
x
y-
x
v
y0
v
x0
g
2
x
2
2
v
x0
Giancoli notes the following. The equation is based on equation 3-6.
We see that as a function of has the form
y
x
y==Ax-B
2
x
yAx-B
2
x
where and are constants for any specific projectile motion. This is the standard equation for a parabola.
A
B
I would like to mention that there should be a constant term
y==Ax-B+
2
x
1
yAx-B+
2
x
c
1
We can rearrange this in decreasing order of the exponent for x by rewriting the equation.
y==A-B+
1
x
2
x
1
0
x
yA-B+
1
x
2
x
c
1
0
x
Rearranging:
Giancoli also notes the following in section *3.7 Projectile Motion Is Parabolic (the star * indicates the section is optional. Giancoli notes * Some Sections of this book, such as this one, may be considered optional at the discretion of the instructor. See the Preface for more details.) on page 64.
considerations
considerations
This equation would be complicated by the following considerations:
geoid
geoid
The gravitational field is not constant. This requires modelling the Earth as a geoid. The best geoid is the Earth Gravitational Model 2020 (EGM2020) released by the National Geospatial-Intelligence Agency (NGA). For more information, see the following links.
I would like to share the introduction:
The National Geospatial-Intelligence Agency [NGA], in conjunction with its U.S. and international partners, has completed its next Earth Gravitational Model (EGM2020), to replace EGM2008. The new ‘Earth Gravitational Model 2020’ [EGM2020] will retain the same harmonic basis and resolution as EGM2008. As such, EGM2020 will be a ellipsoidal harmonic model up to degree (n) and order (m) 2159, but will be released as a spherical harmonic model to degree 2190 and order 2159. EGM2020 has benefited from new data sources and procedures. Updated satellite gravity information from the GOCE and GRACE mission, will better support the lower harmonics, globally. Multiple new acquisitions (terrestrial, airborne and shipborne) of gravimetric data over specific geographical areas (Antarctica, Greenland …), will provide improved global coverage and resolution over the land, as well as for coastal and some ocean areas. Ongoing accumulation of satellite altimetry data as well as improvements in the treatment of this data, will better define the marine gravity field, most notably in polar and near-coastal regions. NGA and partners are evaluating different approaches for optimally combining the new GOCE/GRACE satellite gravity models with the terrestrial data. These include the latest methods employing a full covariance adjustment. NGA is also working to assess systematically the quality of its entire gravimetry database, towards correcting biases and other egregious errors where possible, and generating improved error models that will inform the final combination with the latest satellite gravity models. Outdated data gridding procedures have been replaced with improved approaches. For EGM2020, NGA intends to extract maximum value from the proprietary data that overlaps geographically with unrestricted data, whilst also making sure to respect and honor its proprietary agreements with its data-sharing partners.
Wikipedia notes:
EGM2020 is to be a new release (still not released as of February 2022) with the same structure as EGM2008, but with improved accuracy by incorporating newer data. It was originally planned to be released in April 2020. The precursor version XGM2016 (X stands for experimental) was released in 2016 up to degree and order (d/o) 719. XGM2019e was released in 2020 up to spheroidal d/o 5399 (that corresponds to a spatial resolution of 2′ which is ~4 km) and spherical d/o 5540 with a different spheroidal harmonic construction followed by conversion back into spherical harmonics. XGM2020 was also released recently.
air resistance effect
air resistance effect
Take air resistance into account. This involves atmospheric dynamics. There are several functions that can do this in the Wolfram Language including StandardAtmosphereData which is based on the 1976 United States Standard Atmosphere model.
There are several other models including
ICAO Standard Atmosphere published in 1993 (Mathematica supports the one published in 1976)
weather effects
weather effects
There are also other weather models, including Numerical Weather Prediction models (NWP models). Here are some important ones:
These weather models could be used to model the winds and temperatures that might affect an object.
magnetic effects
magnetic effects
There are also magnetic models, including
Magnetic models could be important if the object in the air is composed of metal.
parabolic motion
parabolic motion
Here is a picture of parabolic motion:
solution
solution
I don't know what the solution to the exercise is. Please comment if you have an idea on how to figure out how high the wall of the fort should be to guarantee that a cannon ball does not go over the wall. The answer should be obtained with one or more differential equations that can be solved with initial conditions to produce an output in terms of L, the distance the cannon is from the fort measured horizontally, H, the height of the hill that the cannon is on, and v, the muzzle velocity of the cannonball. The angle θ that the cannon launches the cannonball might also be significant.
projectile motion
projectile motion
There is lots of helpful information on the Wikipedia page for projectile motion. I would like to analyze the situation where you have to take the curvature of the Earth into account when you have a range greater than 100 km, for example launching an ICBM.
Exercise 8
Exercise 8
I kind of already solved this exercise.
From page 20 and 21.
Solve the IVP and show that the height is given by the familiar formula
Its kind of hard to tell that these are code equivalent.
Exercise 3
Exercise 3
In[]:=
x''(t)=-3sqrt(t),x(1)=1,x'(1)=2
Exercise 1
Exercise 1
Maybe Wolfram Alpha will plot the solution.
In[]:=
x'=t cos t^2
In[]:=
x'=t cos t^2, x(0)=1
Exercise 4
Exercise 4
Find all functions that solve the differential equations:
In[]:=
x'=t E^(-2*t)
In[]:=
sqrt(t) x'=cos sqrt(t)
Exercise 7
Exercise 7
Exercise 11
Exercise 11
A function is defined by the expression
I'm stuck on this one.