Help with solving first order ODE using a simple Fortran code, please

Click For Summary
SUMMARY

This discussion focuses on solving a first-order ordinary differential equation (ODE) using Fortran code. The equation is defined as $$ ds/dt=k_i * \sqrt{v} $$, with variables $$ k_i $$ and $$ v $$ dependent on height $$ h $$, where $$ k_i=\sqrt{χ/h^2} $$ and $$ v= \mu h $$. The user, Joseph, seeks to plot $$ s $$ as a function of $$ h $$ over time but encounters issues with generating multiple plots due to nested loops in the code. The solution involves swapping the loops to allow for a single plot of $$ s $$ against $$ h $$ as time progresses.

PREREQUISITES
  • Understanding of Fortran programming language
  • Familiarity with numerical methods for solving ODEs
  • Knowledge of plotting data in Fortran
  • Basic concepts of differential equations
NEXT STEPS
  • Learn how to implement nested loops in Fortran for time-dependent simulations
  • Research techniques for plotting data using Fortran libraries
  • Explore numerical integration methods for ODEs, such as Euler's method
  • Study the impact of parameter changes on the behavior of differential equations
USEFUL FOR

This discussion is beneficial for Fortran programmers, engineers, and researchers working with numerical simulations of differential equations, particularly those interested in plotting results over time.

joseph2015
Messages
14
Reaction score
2
New user has been reminded to use the Template when posting schoolwork problems
I am trying to solve the following first order ODE using a simple Fortran code :

$$ ds/dt=k_i * \sqrt{v}$$

where both (ki) and (v) are variables depending on (h) as follows

$$ k_i=\sqrt{χ/h^2}$$
$$v= \mu h$$

where (μ) and (χ) are constants. (the arbitrary values of each of them can be seen from the code attached).

Note: h changes from 100 to 1 with -1 reduction

Now after writing the code, I would like to plot (s) as a function of (h). but what I will get is not a single plot, but a series of plots depending on the number of timesteps. the issue has to do with loops, but I am not sure how to get a single plot of (s) against (h) as it grows over a period of time.

Fortran:
program height
    implicit none

    double precision ki, t, v, s, a, x, tend, dt
    integer count,h,nstep
    real mu
    REAL, PARAMETER :: Pi = 3.1415927
    tend   = 100
    dt     = 10
    nstep  = tend/dt

    write(*,*) 'end time for the integration: ', tend
    write(*,*) 'integration time step: ', dt
    write(*,*) 'number of time steps: ', nstep

    open(30, file='height.agr')
 
    x      = 1.99
    s      = 1.0e-2
    mu     = 2.34
    t      = 0.0

    write(30,*) s, real(h*1.8e+4)

    ! START to Integrate.

    write(*,*) 'Start time integration (', nstep, ' steps) ...'

    do count = 1, nstep
        ! update the time
        t = t + dt

!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
       do h  =100,1,-1
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 
          v  =mu*h
          ki =sqrt((x)/h**2)
          a  = ki*sqrt(v)
          s  = s+dt*a

            write(30,*) s, real(h)
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
       end do
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    end do

    write(*,*) 'end'
    close(30)
    100  format(d14.8,1x,d14.8)
end program height
 

Attachments

  • sh.png
    sh.png
    5.9 KB · Views: 356
Last edited:
Physics news on Phys.org
Your write statement is in an inner loop, so you get 100 writes for every time step.
What's worse, s is incremented 100 times for every time step too.
And the outer loop is done only 10 times.

What is it you want to do ?
 
  • Like
Likes   Reactions: berkeman
BvU said:
Your write statement is in an inner loop, so you get 100 writes for every time step.
What's worse, s is incremented 100 times for every time step too.
And the outer loop is done only 10 times.

What is it you want to do ?

Thank you very much for your answer

I am trying to follow the change of s as a function of time and h, then plot s vs h.
for example: if s is a strength of a material, while h is height (from 100 to 1cm), I would like to plot s vs h for 100 years time.

Kind regards,

Joseph
 
My guess is you want to swap the loops and reset t to 0 in the outer loop. That way you get s(t) with h as a parameter.
And I'd do 100 steps in t and 10 steps in h.

And with double $ to delimit the ##\TeX## source you get displayed math:
$$ \frac {ds}{dt}=k_i \sqrt v$$ $$ k_i=\sqrt{\chi/h^2}=\frac {\sqrt \chi }{h } $$
 
BvU said:
My guess is you want to swap the loops and reset t to 0 in the outer loop. That way you get s(t) with h as a parameter.
And I'd do 100 steps in t and 10 steps in h.

And with double $ to delimit the ##\TeX## source you get displayed math:
$$ \frac {ds}{dt}=k_i \sqrt v$$ $$ k_i=\sqrt{\chi/h^2}=\frac {\sqrt \chi }{h } $$

I will try that, thank you very much
 

Similar threads

  • · Replies 4 ·
Replies
4
Views
2K
Replies
2
Views
1K
  • · Replies 4 ·
Replies
4
Views
3K
  • · Replies 5 ·
Replies
5
Views
2K
  • · Replies 7 ·
Replies
7
Views
3K
Replies
18
Views
3K
  • · Replies 2 ·
Replies
2
Views
7K
Replies
7
Views
2K
  • · Replies 12 ·
Replies
12
Views
2K
  • · Replies 5 ·
Replies
5
Views
2K