Proyecto Tronador Aprox

!Este Programa fue realizado por Lillieskold, Thomas Enrique
!Instituto Universitario Aeronautico - 2023
!Docentes: Giovacchini, Juan Pablo; Sacco, Carlos
!Asignatura: Calculo Numerico
!-------------------------------------------------------------------------------
MODULE Datos
REAL(8),PARAMETER:: t_min=0.541
REAL(8),PARAMETER:: t_max=9.d0
REAL(8),PARAMETER:: k=5e-4
REAL(8),PARAMETER:: g=9.80665
REAL(8),PARAMETER:: C_d=0.5
REAL(8),PARAMETER:: rho=1.225
REAL(8),PARAMETER:: s=0.05
REAL(8),PARAMETER:: phi_1=10.d0*sqrt(3.d0)/2.d0
REAL(8),PARAMETER:: phi_2=37.33*sqrt(3.d0)/2.d0
REAL(8),PARAMETER:: phi_3=5.d0
REAL(8),PARAMETER:: phi_4=37.33/2.d0
INTEGER,PARAMETER:: dim =((t_max-t_min)/k) + 25
INTEGER,PARAMETER:: neq=4
END MODULE

Program Euler_explicito_vectorial
USE DATOS
IMPLICIT NONE
REAL(8),DIMENSION(dim,neq)::Numerical_sol

Numerical_sol(1,1)=phi_1
Numerical_sol(1,2)=phi_2
Numerical_sol(1,3)=phi_3
Numerical_sol(1,4)=phi_4

CALL Euler_exp_vec(Numerical_sol)

END PROGRAM Euler_explicito_vectorial

SUBROUTINE Euler_exp_vec(Vn)
USE DATOS
IMPLICIT NONE
REAL(8),DIMENSION(dim,neq)::Vn
REAL(8),DIMENSION(neq)::f
REAL(8) t_i
INTEGER i,j
i=1
t_i=t_min
OPEN(10,FILE="write_outs_k=0.001")

DO WHILE (t_i .le. t_max)

    CALL Fun(Vn(i,:),f,t_i)

    DO j=1,neq,1
        Vn(i+1,j)=Vn(i,j)+k*f(j)
    END DO
    Write(10,*) t_i,Vn(i,1),Vn(i,2),Vn(i,3),Vn(i,4)

    i=i+1
    t_i=t_i+k

END DO
END SUBROUTINE

SUBROUTINE Fun(aprox,f,t)
USE DATOS
IMPLICIT NONE
REAL(8),DIMENSION(neq)::f
REAL(8),DIMENSION(neq)::aprox
REAL(8) t,cos,sin,m,d,thrust

    thrust=5000
    D= 0.5 *C_d *s * rho * (aprox(2)**2+aprox(4)**2)

    cos=aprox(2)/sqrt(aprox(2)**2+aprox(4)**2)

    sin=aprox(4)/sqrt(aprox(2)**2+aprox(4)**2)

    m= (-40.d0/9.d0)*t + 70.d0

    f(1)=aprox(2)

    f(2)= (-(D-thrust)*cos)/m

    f(3)= aprox(3)

    f(4)= (-(m*g) - sin*(D-thrust))/m

END SUBROUTINE