PROGRAM MYPROG
INTEGER NX, NY ! Variables to hold number of intervals in each direction
REAL XSTART, XEND ! Variables defining the interval in x-direction
REAL YSTART, YEND ! Variables defining the inverval in y-direction
REAL DX, DY
!REAL PHI(21,21)
!REAL EQS((20*20),(20*20)+1)
INTEGER COUNTER ! TEMPORARY, USED FOR PRINTING VALUES
REAL, DIMENSION(: , : ), ALLOCATABLE :: PHI
REAL, DIMENSION(: , : ), ALLOCATABLE :: EQS
CALL DEFINEMATRIX(XSTART, XEND, YSTART, YEND, NX, NY) ! Lets the user specify the interval for which the problem will be solved, and the number of intervals for discretisation
ALLOCATE (PHI(NY+1,NX+1))
ALLOCATE (EQS(NY*NX,(NY*NX)+1))
DX=(XEND-XSTART)/NX
DY=(YEND-YSTART)/NY
CALL ENTERBOUNDARYVALUES(PHI,XSTART,DX,NX,YSTART,DY,NY)
!PRINTING OUT TO CONSOLE, TEMPORARY
COUNTER=1
DO WHILE(COUNTER<=NY+1)
WRITE(*,*) PHI((NY+2)-COUNTER,1:NX+1)
COUNTER=COUNTER+1
END DO
!ENDPRINT
CALL CREATEEQUATIONS(PHI, NX, NY, DX, DY, EQS)
END PROGRAM
!---------------------------------!
SUBROUTINE DEFINEMATRIX(XS,XE,YS,YE,NX,NY)
REAL XS, XE, YS, YE
INTEGER NX, NY
WRITE(*,*) 'Choose the interval in x-direction for which you wish to solve the problem. Input format <x_min> <x_max>'
CALL READINTERVAL(XS, XE) ! Prompt user to input x-interval
WRITE(*,*) 'Choose the interval in y-direction for which you wish to solve the problem. Input format <y_min> <y_max>'
CALL READINTERVAL(YS, YE) ! Prompt user to input y-interval
WRITE(*,*) 'Specify number of intervals you want in the x-direction'
READ(*,*) NX
WRITE(*,*) 'Specify number of intervals you want in the y-direction'
READ(*,*) NY
END SUBROUTINE
!----------------------------------!
SUBROUTINE READINTERVAL(START, END)
REAL START, END
read(*,*) START, END ! Assignes the values to the corresponding variables
END SUBROUTINE
!----------------------------------!
REAL FUNCTION G(X,Y)
REAL X,Y
G=1
RETURN
END FUNCTION
!----------------------------------!
REAL FUNCTION PHIBOUNDARY(X,Y)
REAL X, Y
PHIBOUNDARY=(X**2+Y**2)/4.0
RETURN
END FUNCTION
!----------------------------------!
SUBROUTINE ENTERBOUNDARYVALUES(PHI, XSTART, DX, NX, YSTART, DY, NY)
REAL PHI(NY+1,NX+1)
REAL XSTART, DX, YSTART, DY
INTEGER NX, NY
INTEGER COL
INTEGER ROW
DO WHILE(COL<=NX+1)
PHI(1,COL)=PHIBOUNDARY(XSTART+(COL-1)*DX,YSTART)
PHI(NY+1,COL)=PHIBOUNDARY(XSTART+(COL-1)*DX,YSTART+NY*DY)
COL=COL+1
END DO
ROW=2
DO WHILE(ROW<=NY)!DENNE GÅR FRA 2 TIL OG MED NY, DA 1 OG NY+1 ER FYLT UT AV LØKKEN OVER
PHI(ROW,1)=PHIBOUNDARY(XSTART,YSTART+(ROW-1)*DY)
PHI(ROW,NX+1)=PHIBOUNDARY(XSTART+NX*DX,YSTART+(ROW-1)*DY)
ROW=ROW+1
END DO
END SUBROUTINE
SUBROUTINE CREATEEQUATIONS(PHI, NX, NY, DX, DY, EQS)
REAL PHI(NY+1,NX+1)
INTEGER NX, NY
REAL DX, DY
REAL EQS(NX*NY,(NX*NY)+1)
REAL EQ((NX*NY)+1)
INTEGER COUNTX, COUNTY
INTEGER ROW, COL
REAL COEFFDY, COEFFDX, COEFFG
COEFFDY=-1.0*(DX**2)/(2.0*((DX**2)+(DY**2)))
COEFFDX=-1.0*(DY**2)/(2.0*((DX**2)+(DY**2)))
COEFFG=(DX**2)*(DY**2)/(2.0*((DX**2)+(DY**2)))
COUNTY=2
EQS=0.0
EQ=0.0
DO WHILE(COUNTY<=NY)
COUNTX=2
DO WHILE(COUNTX<=NX)
!Antar at vi befinner oss på en av de indre indre nodene, (sjekk med if-setninger)
ROW=COUNTY
COL=COUNTX
!WRITE(*,*) '(',COL,',',ROW,')'
!WRITE(*,*) INT(GET1DPOS(COL, ROW, NX)),' er posisjonen nå'
!WRITE(*,*) TEST(COL, ROW, NX),' er posisjonen nå'
!WRITE(*,*) 'Oppdatert 2'
EQS(COUNTY, INT(GET1DPOS(COL, ROW, NX)) )=1.0
WRITE(*,*) '1'
EQS(COUNTY, INT(GET1DPOS(COL+1,ROW,NX)) )=COEFFDX
WRITE(*,*) '2'
EQS(COUNTY, INT(GET1DPOS(COL-1,ROW,NX)) )=COEFFDX
WRITE(*,*) '3'
EQS(COUNTY, INT(GET1DPOS(COL,ROW+1,NX)) )=COEFFDY
WRITE(*,*) '4'
EQS(COUNTY, INT(GET1DPOS(COL,ROW-1,NX)) )=COEFFDY
EQS(COUNTY, INT(GET1DPOS(2,NX-1,NX)) )=COEFFG
WRITE(*,*) '6'
EQS(COUNTY,(INT(2.2)))=0.55
WRITE(*,*) EQS(COUNTY,:)
!WRITE(*,*) EQ(1:(NX*NY)+1)
!WRITE(*,*), COUNTX,'<=', NX, ' X'
COUNTX=COUNTX+1
!WRITE(*,*) HI(COL,ROW,NX), 'hei'
!WRITE(*,*) 'HEI
END DO
!WRITE(*,*), COUNTY,'<=', NY, ' Y'
COUNTY=COUNTY+1
END DO
!WRITE(*,*) EQS
!TEMP
!COUNTX=1
!COUNTY=1
!WRITE(*,*) 'Dette er en test av lignignsmatrisa'
!DO WHILE(COUNTY<=NY+1)
!WRITE(*,*) EQS(COUNTY,1:NX+2)
!COUNTY=COUNTY+1
!END DO
!ENDTEMP
END SUBROUTINE
REAL FUNCTION GET1DPOS(COL, ROW, NX)
INTEGER COL, ROW, NX
!WRITE(*,*) ' Jeg har kommet hit'
GET1DPOS=(COL-1)+((ROW-2)*(NX-1))
!WRITE(*,*) COL, ROW, NX, GET1DPOS, 'hei'
!WRITE(*,*) GET1DPOS
RETURN
END FUNCTION