! Nicholas Mai
! PHYS 401 - Computational Physics
! Homework 1 - Population Growth

! Linear ODE Version
! This program numerically computes the solution to a first-order, linear
! ordinary differential equation. The method used is the basic Euler explicit
! scheme.

PROGRAM lin_popgrowth
    IMPLICIT NONE

    ! Declaration of arrays
    REAL, DIMENSION(:), ALLOCATABLE:: N,NA,R,RA,t
    REAL:: gamma_,dt,N0
    INTEGER::t_f, num_it, i
    CHARACTER (LEN=20):: n_out
    CHARACTER (LEN=20):: na_out
    CHARACTER (LEN=20):: r_out
    CHARACTER (LEN=20):: ra_out

    ! Ask for input for variables
    WRITE(*,*) 'Enter the final time (t_f): '
        READ(*,*) t_f
    WRITE(*,*) 'Enter the time step (dt): '
        READ(*,*) dt
    WRITE(*,*) 'Enter the gamma factor (gamma): '
        READ(*,*) gamma_
    WRITE(*,*) 'Enter the computed N(t) file name (n_out): '
        READ(*,*) n_out
    WRITE(*,*) 'Enter the computed R(t) file name (r_out): '
        READ(*,*) r_out
    WRITE(*,*) 'Enter the analytical N(t) file name (na_out): '
        READ(*,*) na_out
    WRITE(*,*) 'Enter the analytical N(t) file name (ra_out): '
        READ(*,*) ra_out

    ! Input / Output parameters
    OPEN(7,FILE=n_out)
    OPEN(8,FILE=r_out)
    OPEN(9,FILE=na_out)
    OPEN(10,FILE=ra_out)

    ! Number of iterates
    num_it = INT(t_f/dt)
    WRITE(*,*) 'Number of iterations: '
    WRITE(*,*) num_it

    ! Allocate and initialize the main variables
    ALLOCATE(N(0:(num_it-1)))
    ALLOCATE(R(0:(num_it-1)))
    ALLOCATE(t(0:(num_it-1)))
    WRITE (*,*) 'Enter the initial value: '
        READ(*,*) N0
    N(0) = N0
    R(0) = N0 * gamma_
    t(0) = 0
    NA(0) = N0
    RA(0) = N0 * gamma_

    ! Calculate the solution via the explicit Euler's method
    DO i = 0, num_it, 1
        N(i+1) = N(i) +  N(i) * gamma_ * dt;
        R(i+1) = N(i) * gamma_;
        NA(i) = N0 * exp(gamma_ * t(i));
        RA(i) = gamma_ * NA(i);
        t(i+1) = t(i) + dt;
        WRITE(7,*)t(i),N(i+1)
        WRITE(8,*)t(i),R(i+1)
        WRITE(9,*)t(i),NA(i)
        WRITE(10,*)t(i),RA(i)
    END DO

    ! Closing all input / output connections
    CLOSE(7)
    CLOSE(8)
    CLOSE(9)
    CLOSE(10)
STOP
END