# SymPy to Solve Initial Value Condition (First Order Differential Equation)

**URL:** <https://discourse.julialang.org/t/sympy-to-solve-initial-value-condition-first-order-differential-equation/97203>\
**Category:** General Usage\
**Tags:** question\
**Created:** [April 7, 2023, 4:08am UTC](https://discourse.julialang.org/t/sympy-to-solve-initial-value-condition-first-order-differential-equation/97203 "2023-04-07T04:08:29Z")\
**Posts on this page:** 3\
**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 7, 2023, 4:08am UTC](https://discourse.julialang.org/t/sympy-to-solve-initial-value-condition-first-order-differential-equation/97203/1 "2023-04-07T04:08:29Z")

</div>

Hi all,

I want to check is my code correct to compute IVP of

\frac{dx}{dt} = \sqrt{2g} \frac{R}{\sqrt{R+x}}

the initial condition is x(0)=0

```julia
using SymPy
@syms t, g, R
@syms x()

# Computing IVP with SymPy
x(t).diff(t) 
diffeq = Eq(x(t).diff(t), sqrt(2g)*R*(1/sqrt(R+x(t))) ); string(diffeq)

# To solve the ODE, pass it and the function to solve for to dsolve.
xt = dsolve(diffeq, x(t))
diffeq3 = xt[3](Dict(R => 4000, g => 78545 ))
# To solve the Initial Value Problem, initial condition is 0
ivp = dsolve(diffeq3, x(t), ics = Dict(x(0) => 0 ))

```

there are 3 solutions at xt and when I try to compute IVP it becomes error:

**PyError ($(Expr(:escape, :(ccall(#= /home/browni/.julia/packages/PyCall/ygXW2/src/pyfncall.jl:43 =# @pysym(:PyObject\_Call), PyPtr, (PyPtr, PyPtr, PyPtr), o, pyargsptr, kw))))) \<class ‘ValueError’\>**  
**ValueError(“Couldn’t solve for initial conditions”)**  
\*\* File “/home/browni/.julia/conda/3/envs/lasthrim\_env/lib/python3.10/site-packages/sympy/solvers/ode/ode.py”, line 640, in dsolve\*\*  
\*\* return \_helper\_simplify(eq, hint, hints, simplify, ics=ics)\*\*  
\*\* File “/home/browni/.julia/conda/3/envs/lasthrim\_env/lib/python3.10/site-packages/sympy/solvers/ode/ode.py”, line 695, in \_helper\_simplify\*\*  
\*\* solved\_constants = solve\_ics([rv], [r[‘func’]], cons(rv), ics)\*\*  
\*\* File “/home/browni/.julia/conda/3/envs/lasthrim\_env/lib/python3.10/site-packages/sympy/solvers/ode/ode.py”, line 802, in solve\_ics\*\*  
\*\* raise ValueError(“Couldn’t solve for initial conditions”)\*\*

**Stacktrace:**

---

<div class="post-metadata">

**Author:** ![j\_verzani](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/j_verzani/32/8551_2.png) [@j\_verzani](https://discourse.julialang.org/u/j_verzani)\
**Post date:** [April 7, 2023, 12:24pm UTC](https://discourse.julialang.org/t/sympy-to-solve-initial-value-condition-first-order-differential-equation/97203/2 "2023-04-07T12:24:06Z")

</div>

This isn’t quite the way to set this up. You are calling `dsolve` twice, not once. However, here SymPy also gives up. You can do the work manually, though it isn’t so pretty (sorry, I switched some variable names):

```julia
julia> eq = diff(u(t), t) - √(2g) * R / √(R + u(t))

julia> sols = dsolve(eq);

julia> sol = rhs(sols[3]);

julia> c1 = free_symbols(sol)[1]
C₁

julia> C1 = solve(sol(t=>0) ~ 0, c1)[2];

julia> sol(c1 => C1)
...

```

---

<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 8, 2023, 5:31am UTC](https://discourse.julialang.org/t/sympy-to-solve-initial-value-condition-first-order-differential-equation/97203/3 "2023-04-08T05:31:57Z")

</div>

Yes thanks @j_verzani , I am starting to learn how to use SymPy to solve DE now.

this is working code:

```julia
using SymPy
@syms t, g, R
@syms x()

# Computing IVP with SymPy
diffeq = diff(x(t), t) - √(2g) * R / √(R + x(t))

# To solve the ODE, pass it and the function to solve for to dsolve.
xt = dsolve(diffeq)

# To solve the Initial Value Problem, initial condition is 0
sol_xt = rhs(xt[3])
c1 = free_symbols(sol_xt)[1]
t1 = free_symbols(sol_xt)[4]

C1 = solve(sol_xt(t=>0) ~ 0, c1)[2]
T1 = solve(sol_xt(c1=>C1) + 4000 ~ 240000, t1)[2]

ivp = sol_xt(c1=>C1)
#ivp = sol_xt(c1=>C1,R => 4000, g => 78545)

println("Modeling with First Order Equations with Julia")
println("A Rocket is Launched straight up from the surface of the earth")
println("Neglect air resistance")
println("The initial condition: x(0) = 0")
println("g = GM/(R^2)")
println("g = 78,545 mi/hr^2")

println("")
println("The time required for the rocket to go 240,000 miles (in hours): ")
println(T1(R => 4000, g => 78545).evalf())

println("")
println("The initial value problem solution is:")
println("x(t) = ")
ivp

```