# Checking the Solutions for First Order Differential Equations

**URL:** <https://discourse.julialang.org/t/checking-the-solutions-for-first-order-differential-equations/97304>\
**Category:** General Usage\
**Tags:** question\
**Created:** [April 10, 2023, 1:39pm UTC](https://discourse.julialang.org/t/checking-the-solutions-for-first-order-differential-equations/97304 "2023-04-10T13:39:47Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [April 10, 2023, 1:39pm UTC](https://discourse.julialang.org/t/checking-the-solutions-for-first-order-differential-equations/97304/1 "2023-04-10T13:39:47Z")

</div>

Hi all,

I just want to make sure, whether the solution of x(t) and y(t) are correct in t term.

**If the solutions are incorrect, it is `SymPy` fault.**

Because if I put t\>0 I get x(t) \< 0 and I think it should not be possible. x(t) is the horizontal coordinate while y(t) is the vertical coordinate of a projectile motion.

the code:

```julia
using SymPy, Plots
@syms t, g, r, A, u, h
@syms v()
@syms w()
@syms x()
@syms y()

# Computing IVP with SymPy
diffeqv = Eq(v(t).diff(t), -r*v(t)); string(diffeqv)
diffeqw = Eq(w(t).diff(t), -g -r*w(t)); string(diffeqw)

# To solve the ODE, pass it and the function to solve with dsolve.
vt = dsolve(diffeqv, v(t))
wt = dsolve(diffeqw, w(t))

sol_vt = rhs(vt)
sol_wt = rhs(wt)

# To solve the Initial Value Problem, initial condition v(0) = u cos (A)
ivp_vt = dsolve(diffeqv, v(t), ics = Dict(v(0) => u*cos(A) ))
ivp_wt = dsolve(diffeqw, w(t), ics = Dict(w(0) => u*sin(A) ))

# Manual way
# Find the constant C₁
c1 = free_symbols(vt)[1]
# To solve the Initial Value Problem, initial condition v(0) = u cos (A)
C1 = solve(sol_vt(t=>0) ~ u*cos(A), c1)
#ivp_vt = sol_vt(c1=>C1[1])

# Part (b)
diffeqx = Eq(x(t).diff(t), rhs(ivp_vt));
diffeqy = Eq(y(t).diff(t), rhs(ivp_wt));

# To solve the ODE, pass it and the function to solve with dsolve.
xt = dsolve(diffeqx, x(t))
yt = dsolve(diffeqy, y(t))

sol_xt = rhs(xt)
sol_yt = rhs(yt)
ivp_xt = dsolve(diffeqx, x(t), ics = Dict(x(0) => 0 ))
ivp_yt = dsolve(diffeqy, y(t), ics = Dict(y(0) => h ))

println("Modeling with First Order Equations with Julia")
println("A baseball is thrown with the initial speed of u and initial angle of elevation is A")
println("v(t) is the horizontal component of the velocity of the thrown baseball")
println("w(t) is the vertical component of the velocity of the thrown baseball")
println("x(t) is the horizontal coordinate of the thrown baseball at time t")
println("y(t) is the vertical coordinate of the thrown baseball at time t")

sympy.pretty_print(ivp_vt)
sympy.pretty_print(ivp_wt)
println("")
println("The horizontal and vertical coordinates : ")
sympy.pretty_print(ivp_xt)
sympy.pretty_print(ivp_yt)

```

![Capture d’écran_2023-04-10_20-38-20](https://global.discourse-cdn.com/julialang/original/3X/5/5/5504bce46db9871f4a46fcd8414f6321da542d7a.png)

To check via PrettyTables:

```julia
using SymPy, Plots, LaTeXStrings, Plots.PlotMeasures, PrettyTables

gr()
u = 125
r = 1/5
h = 3
A = 15  
g = 32

#int = range(0, stop=10, length=300)
int = -3π:1:0.2π
x(t) = (u/r)*cos(A) - (u/r)*exp(-r*t)*cos(A)
# 625*cos(A)-625*exp(-r*t)*cos(A)
y(t) = -(g*t/r) + (g + h*r^2 + r*u*sin(A))/(r^2) - (u*sin(A)/r + (g/r^2))*exp(-r*t)

println(pretty_table(String, 
 	 vcat([[string(t) x(t) y(t)] for t in -4:0.5:4]...); 
  	header=["t", "x(t)", "y(t)"]))

```

![Capture d’écran_2023-04-10_20-39-14](https://global.discourse-cdn.com/julialang/original/3X/7/2/728b58c6237c251205f728e843640d5807778025.png)

The solutions of x(t) and y(t) from the book’ solution manual:

![Capture d’écran_2023-04-10_20-41-55](https://global.discourse-cdn.com/julialang/original/3X/5/d/5d40dbb1b0bb5f33fc78adaa08354682336fc96c.png)

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [April 10, 2023, 2:37pm UTC](https://discourse.julialang.org/t/checking-the-solutions-for-first-order-differential-equations/97304/2 "2023-04-10T14:37:57Z")

</div>

I get the desired plot by choosing A=20,45,70

```julia
using SymPy, Plots, LaTeXStrings, Plots.PlotMeasures

gr()
u = 125
r = 1/5
h = 3
A1 = 20  
A2 = 45  
A3 = 70  
g = 32

#int = range(0, stop=10, length=300)
int = 0:0.1:2π
x1(t) = (u/r)*cos(A1) - (u/r)*exp(-r*t)*cos(A1)
x2(t) = (u/r)*cos(A2) - (u/r)*exp(-r*t)*cos(A2)
x3(t) = (u/r)*cos(A3) - (u/r)*exp(-r*t)*cos(A3)

y1(t) = -(g*t/r) + (g + h*r^2 + r*u*sin(A1))/(r^2) - (u*sin(A1)/r + (g/r^2))*exp(-r*t)
y2(t) = -(g*t/r) + (g + h*r^2 + r*u*sin(A2))/(r^2) - (u*sin(A2)/r + (g/r^2))*exp(-r*t)
y3(t) = -(g*t/r) + (g + h*r^2 + r*u*sin(A3))/(r^2) - (u*sin(A3)/r + (g/r^2))*exp(-r*t)

# note the dot, which means compute x for all values in the int range
plot(x1.(int),y1.(int), ylims=(0,200),
	label=L"A=20",
	legend=:topright)
plot!(x2.(int),y2.(int), ylims=(0,200),
	label=L"A=45",
	legend=:topright)
plot!(x3.(int),y3.(int), ylims=(0,200),
	label=L"A=70", left_margin=10mm,
	legend=:topright)
xlabel!(L"x(t) = \frac{u}{r} \cos \ A (1 - e^{-rt}) ")
ylabel!(L"y(t) = - \frac{gt}{r} + \frac{g + ur \ \sin \ A + hr^{2}}{r^{2}} - (\frac{u}{r}")
#=
@syms t
solve(diff(y(t),t)) # to find extreme point

=#

```

---

<div class="post-metadata">

**Author:** ![JM\_Beckers](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jm_beckers/32/22482_2.png) [@JM\_Beckers](https://discourse.julialang.org/u/JM_Beckers)\
**Post date:** [April 11, 2023, 5:58pm UTC](https://discourse.julialang.org/t/checking-the-solutions-for-first-order-differential-equations/97304/3 "2023-04-11T17:58:57Z")

</div>

I think you should realize that cos() etc are working using radians, not degrees.

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [April 17, 2023, 5:46am UTC](https://discourse.julialang.org/t/checking-the-solutions-for-first-order-differential-equations/97304/4 "2023-04-17T05:46:21Z")

</div>

Yes, this reminds me how I should use pi instead of angles… Thanks.