ImplicitTemp

!*****************************************************************************
!
!! This program is to solve one dimensional heat transfer in a rod with constant tempreture boundary conditions by implicit method
!

program oneDTempIm
implicit none
real, allocatable :: T(:), B(:), Lx(:)
 
real :: Dx, R, Dt, Alfa, M, TRight, TLeft, initial, time, u, stability
integer :: i, j, n, nx

!****************************************************

!**************** Introducing Parameters**************
do
print *, "please enter the numbers of Dx "
read *,   Dx
print *, "please enter the numbers of Dt and time of solution "
read *,   Dt, time
print *, "please enter the numbers of T(Left), T(Right)"
read *,   TLeft, TRight
print *, "please enter the numbers of alfa and initial condition"
read *,   alfa, initial
stability=(alfa*Dt)/(Dx**2)
if(stability<0.5) exit
end do

nx = ((1.0/Dx)+1)
n = nx

allocate(T(n), B(n), Lx(n))


R = (Alfa * Dt)/(Dx**2)
M = (1-(2*R))
!======================End===============================

    B(1) = TLeft
    B(n) = TRight
    
    do i= 2, n-1
    B(i) = initial
    end do

    
    u = 0

    do while ( u <= time )


    T(1) = TLeft
    T(n) = TRight
    
    do i = 2, n-1
    T(i) = R*B(i+1) + M*B(i) + R*B(i-1)
    end do


    do i = 1 , n
    B(i) = T(i)
    end do

   u = u + Dt

    end do  

!***********showing Temperatures*************
OPEN(UNIT=15,FILE='Temperatures.txt')
do i = 1, n
write(15,*) i, T(i) 
end do

close(15)

!======================End===============================

!*********showing temperature contour in tecplot**********

do i=1, n
Lx(i)= (i-1)*Dx
end do

OPEN(UNIT=20,FILE='Temperature contour.PLT')
WRITE(20,*)'VARIABLES="X","T"'
WRITE(20,*)'ZONE I=',nx,'F=POINT' 
do i=1,n 
WRITE(20,*)Lx(i), T(i)
END DO
close(20)
!======================End===============================

end program oneDTempIm