Fortran Frustrating fortran .exe problem

[Rahman_M_M]

HI,
I am using fortran 90 for simulation. I write a subroutine (CALC_PORO_SOL_BANDED )where i called a linear solver LSLRB.
I called this CALC_PORO_SOL_BANDED subroutine from another two subroutine in a do loop. It works for 1 subroutine without any error but from another subroutine it works 15-20 times and after that i get an error message "Program Name.EXE HAS ENCOUNTERDED A PROBLEM AND NEEDS TO CLOSE.WE ARE SORRY FOR THE INCONVENIENCE", and program terminates.

Anyone can give me any idea? Please help me out of this problem.Rahman M M

 

[mathmate]

Without data and without code, it is hard to imagine what needs to be fixed.
Do your subroutines need a lot of memory? Do they recurse?

 

[Rahman_M_M]

Thanks mathmate. Yes I have used a lot of data. Actually that subroutine i used make a big matrix (5103x5103), and make it banded as (156x5103), and call for linear solver. herte is that code.

SUBROUTINE CALC_PORO_SOL_BANDED (nN,refNode,nE,npE,ELEM,BAND_WIDTH,nN_W_B,Wellbore_Boundary_Node,nN_O_B,Outer_Boundary_Node, &
& NEW_nN_X,X_axis_Node,nN_Y,Y_axis_Node,Permiab,Mue,d_T,CT,Phi,ALPHA,D,Curr_Pw,Prev_Pw,S_DISP_LOAD,S_PRESS_LOAD,S_UxSol,S_UySol,S_PSol)

INTEGER :: i,nN,nE,npE,nN_O_B,nN_W_B,NEW_nN_X,nN_Y,BAND_WIDTH,NUCA,NLCA,ROW,COL
REAL ::ALPHA,Mue,d_T,CT,Phi,Curr_Pw,Prev_Pw,MAX_BANDED_PORO_MASS_MAT

INTEGER :: ELEM(nE,npE),X_axis_Node(NEW_nN_X),Y_axis_Node(nN_Y),Wellbore_Boundary_Node(nN_W_B),&
& Outer_Boundary_Node(nN_O_B)

REAL:: refNode(nN,2),S_UxSol(nN),S_UySol(nN),S_PSol(nN),Permiab(2),D(3,3), &
& S_PRESS_LOAD(nN),S_DISP_LOAD(nN*2)

REAL,DIMENSION(:),ALLOCATABLE ::M_Load,TotSol
REAL,DIMENSION(:,:),ALLOCATABLE::BANDED_PORO_MASS_MAT
ALLOCATE(TotSol(nN*3),M_Load(nN*3),BANDED_PORO_MASS_MAT(2*BAND_WIDTH+1,nN*3))

!INITIALIZED THE DISPLACEMENT
!TotSol=0.0;S_UxSol=0.0;S_UySol=0.0;S_PSol=0.0
NUCA=BAND_WIDTH
NLCA=BAND_WIDTH
ROW=0;COL=0
M_Load=0.0
DO i=1,nN
M_Load(i*3-2)=S_DISP_LOAD(i*2-1)
M_Load(i*3-1)=S_DISP_LOAD(i*2)
M_Load(i*3)=S_PRESS_LOAD(i)
END DO
!THIS METHOD WILL CALCULATE THE MASS MATRIX (156x5013)!
CALL CAL_BANDED_GLOBAL_MASS_MATRIX(nN,refNode,nE,npE,ELEM,BAND_WIDTH,Permiab,Mue,d_T,CT,Phi,Alpha,D,BANDED_PORO_MASS_MAT)


MAX_BANDED_PORO_MASS_MAT=MAXVAL(BANDED_PORO_MASS_MAT)
MAX_BANDED_PORO_MASS_MAT=MAX_BANDED_PORO_MASS_MAT*100000.0



!################################$FOR KNOWN PRESSURE$####################################$# DO i=1,nN_W_B ROW=Wellbore_Boundary_Node(i)*3 COL=ROW M_Load(ROW)=(Curr_Pw-Prev_Pw)*MAX_BANDED_PORO_MASS_MAT BANDED_PORO_MASS_MAT(ROW+NUCA+1-COL,COL)=BANDED_PORO_MASS_MAT(ROW+NUCA+1-COL,COL)+MAX_BANDED_PORO_MASS_MAT END DO !$##############################$FOR KNOWN DISPLACEMENT$#######################################
!FOR X AND Y AXIS NODE EXCEPT FREE NODE
DO i=1,NEW_nN_X !FOR X AXIS NODE EXCEPT FREE NODE FOR FRACTURE Y DISPLACEMENT IS ZERO
ROW=X_axis_Node(i)*3-1
COL=ROW
M_Load(ROW)=0.0 !0.0*MAX_BANDED_PORO_MASS_MAT
BANDED_PORO_MASS_MAT(ROW+NUCA+1-COL,COL)=BANDED_PORO_MASS_MAT(ROW+NUCA+1-COL,COL)+MAX_BANDED_PORO_MASS_MAT
END DO
DO i=1,nN_Y !FOR X AXIS NODE X DISPLACEMENT IS ZERO
ROW=Y_axis_Node(i)*3-2
COL=ROW
M_Load(ROW)=0.0 !0.0*MAX_BANDED_PORO_MASS_MAT
BANDED_PORO_MASS_MAT(ROW+NUCA+1-COL,COL)=BANDED_PORO_MASS_MAT(ROW+NUCA+1-COL,COL)+MAX_BANDED_PORO_MASS_MAT
END DO


CALL LSLRB (nN*3, BANDED_PORO_MASS_MAT, NLCA+NUCA+1, NLCA, NUCA,M_Load , 1, TotSol)
!REARRANGE THE PRESSURE AND DISPLACEMENT
Do i=1,nN
S_UxSol(i)=TotSol(i*3-2)
S_UySol(i)=TotSol(i*3-1)
S_PSol(i)=TotSol(i*3)
End Do
DEALLOCATE(BANDED_PORO_MASS_MAT)
DEALLOCATE(M_Load)
DEALLOCATE(TotSol)

END SUBROUTINE CALC_PORO_SOL_BANDED

This subroutine give me the solution.I call this subroutine from another two module those are called from main.

MODULE PORO_DISPL_WITH_GAUSS_PRESS

USE PORO_BNDRY_LOAD
USE GLOBAL_MASS_MATRIX_BANDED
USE PORO_FRAC_LOAD
USE BASIC_CALC
USE MSIMSL
USE PORTLIB


IMPLICIT NONE
PUBLIC CALC_DISPLMNT_WITH_GAUSS_PRESS
PRIVATE CALC_PORO_SOL_BANDED_2
CONTAINS

!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
SUBROUTINE CALC_DISPLMNT_WITH_GAUSS_PRESS(nN,refNode,nN_W_B,Wellbore_Boundary_Node,n_WB_E, &
& N_TOTAL_FRAC_ELEM,NO_FREE_ELEM,nN_X,FRAC_X_AXIS_NODE,Curr_Pw,Prev_S_Prop_Pw,STEP, &
& Sigma_XX_INIT,Sigma_YY_INIT,REF_LEL_FRAC_ELEM,M2_Matr,Prev_P_P,nE,npE,ELEM,BAND_WIDTH, &
& nN_O_B,Outer_Boundary_Node,NO_New_KN_X_Axis_Node,New_KN_X_Axis_Node,nN_Y,Y_axis_Node, &
& Permiab,Mue,del_T_g,CT,Phi,ALPHA,D,PORO_DISPL_BNDY_LOAD,del_Ux,del_Uy,del_P)


INTEGER :: STEP,nN,nE,npE,BAND_WIDTH,NO_NEW_KN_X_AXIS_NODE,nN_X,nN_Y,NN_O_B,nN_W_B,N_WB_E,N_TOTAL_FRAC_ELEM,NO_FREE_ELEM

INTEGER :: New_KN_X_Axis_Node(NO_NEW_KN_X_AXIS_NODE),FRAC_X_AXIS_NODE(nN_X),Y_axis_Node(nN_Y),OUTER_BOUNDARY_NODE(NN_O_B), &
& WELLBORE_BOUNDARY_NODE(nN_W_B),REF_LEL_FRAC_ELEM(N_TOTAL_FRAC_ELEM,4),ELEM(nE,npE)

REAL :: ALPHA,Prev_S_Prop_Pw,CURR_PW,del_T_g,Permiab(2),Mue,CT,Phi

REAL :: refNode(nN,2),FRACTURED_DISP_LOAD(nN*2),PRESS_LOAD(nN),del_Ux(nN),del_Uy(nN),del_P(nN),Sigma_XX_INIT(nN),Sigma_YY_INIT(nN),M2_Matr(nN,nN),&
& DISP_LOAD(2*nN),PORO_DISPL_BNDY_LOAD(2*nN),Prev_P_P(nN),D(3,3)

!THIS FUNCTION WILL CALCULATE THE DISPL LOAD ALONG BOUNDARY
PORO_DISPL_BNDY_LOAD=0.0
CALL CALC_PORO_ELAS_DISPL_WELLBORE_BNDY_LOAD(nN,refNode,nN_W_B,Wellbore_Boundary_Node,n_WB_E,Prev_S_Prop_Pw,Curr_Pw,PORO_DISPL_BNDY_LOAD)

!THIS FUNCTION WILL CALCULATE THE DISPL LOAD ALONG FRACTURE
FRACTURED_DISP_LOAD=0.0

CALL CALC_PORO_FRACTURED_DISPL_LOAD_WITH_GAUSS_BHP(nN,refNode,N_TOTAL_FRAC_ELEM,NO_FREE_ELEM,nN_X,FRAC_X_AXIS_NODE, &
& Curr_Pw,Prev_S_Prop_Pw,STEP,Sigma_XX_INIT,Sigma_YY_INIT,REF_LEL_FRAC_ELEM,FRACTURED_DISP_LOAD)

!THIS FUNCTION WILL CALCULATE THE PRESSURE LOAD IN ALL NODE
PRESS_LOAD=(-1.0)*del_T_g*MATMUL(M2_Matr,Prev_P_P)
DISP_LOAD=PORO_DISPL_BNDY_LOAD+FRACTURED_DISP_LOAD

PRINT*,"INSIDE THE GAUSS MODEL BEFORE CALL"

CALL CALC_PORO_SOL_BANDED (nN,refNode,nE,npE,ELEM,BAND_WIDTH,nN_W_B,Wellbore_Boundary_Node,nN_O_B,Outer_Boundary_Node, &
& NO_New_KN_X_Axis_Node,New_KN_X_Axis_Node,nN_Y,Y_axis_Node,Permiab,Mue,del_T_g,CT,Phi,ALPHA,D,Curr_Pw,Prev_S_Prop_Pw,DISP_LOAD,PRESS_LOAD,del_Ux,del_Uy,del_P)

PRINT*,"INSIDE THE GAUSS MODEL AFTER CALL ......."
END SUBROUTINE CALC_DISPLMNT_WITH_GAUSS_PRESS

END MODULE PORO_DISPL_WITH_GAUSS_PRESS



Here is the portion of the main where i called that subroutine in a do loop.


*********************************Main*************************
ITERATIVE: DO WHILE(TOLLARENCE.GE.TOLLARENCE_ACC)
del_T_g=0.1
Tol_Pw_Min_H=1.0-(MinHStr/Curr_Pw)
IF(Tol_Pw_Min_H.LE.0.005) THEN
PRINT*,"FRACTURE PRESSURE IS BELLOW THE MINIMUM IN SITU HORIZONTAL STRESS"
PRINT*," PROPAGATION STOP"
EXIT FRC_ELM
END IF
PRINT*,"BEFORE GAUSS MODEL"
!Curr_Pw=Pw

!THIS SUBROUTINE WILL CALCULATE ALL LOAD AND SOLVED THE DISPLACEMENT AND PRESSURE WITH GAUSS BHP
CALL CALC_DISPLMNT_WITH_GAUSS_PRESS(nN,refNode,nN_W_B,Wellbore_Boundary_Node,n_WB_E, &
& N_TOTAL_FRAC_ELEM,NO_FREE_ELEM,nN_X,FRAC_X_AXIS_NODE,Curr_Pw,Prev_S_Prop_Pw,STEP, &
& Sigma_XX_INIT,Sigma_YY_INIT,REF_LEL_FRAC_ELEM,M2_Matr,Prev_P_P,nE,npE,ELEM,BAND_WIDTH, &
& nN_O_B,Outer_Boundary_Node,NO_New_KN_X_Axis_Node,New_KN_X_Axis_Node,nN_Y,Y_axis_Node, &
& Permiab,Mue,del_T_g,CT,Phi,ALPHA,D,PORO_DISPL_BNDY_LOAD,del_Ux,del_Uy,del_P)

PRINT*,"AFTER GAUSS MODEL"
PRINT*,"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"

UX_WITH_G_PRESS=del_Ux+UX_INIT
UY_WITH_G_PRESS=del_Uy+UY_INIT

OPEN(1111, FILE ="PORO_PRESSURE_WITH_FLATP.txt")
WRITE(1111,1)del_P*Pres_SCfactor/PSI_TO_PA
CLOSE(1111)

OPEN(1010, FILE ="PORO_XDIDP_WITH_FLATP.txt")
WRITE(1010,1)UX_WITH_G_PRESS
CLOSE(1010)


!PRINT*,"TEST"
!PAUSE

!WIDTH_WITH_G_PRESS=0.0
DO k=1,nN_X
WIDTH_WITH_G_PRESS(k)=Uy_WITH_G_PRESS(FRAC_X_AXIS_NODE(k))
END DO

!refNode_Upd=0.0
refNode_Upd(:,1)=refNode(:,1)-Ux_WITH_G_PRESS
refNode_Upd(:,2)=refNode(:,2)-Uy_WITH_G_PRESS


!CALL FLUID FLOW INSIDE THE FRACTURE AND DISPLACEMENT CONVERGENCE FOR FRACTURE PRESSURE
PRINT*,"BEFORE EXACT MODEL"
CALL CALC_DISPL_CONV_AND_FRAC_FLUID_FLOW_PORO(nN,refNode,refNode_Upd,NO_New_KN_X_Axis_Node,New_KN_X_Axis_Node, &
& nN_Y,Y_axis_Node,nN_W_B,Wellbore_Boundary_Node,nN_O_B,Outer_Boundary_Node,nE,npE,ELEM,BAND_WIDTH,nN_X, &
& FRAC_X_AXIS_NODE,N_TOTAL_FRAC_ELEM,NO_FREE_ELEM,LEAKOFF_COFF,Meu_FR_F,QRATE,FRAC_WIDTH_PREV,Uy_WITH_G_PRESS, &
& del_T,PREV_NODAL_PRESS_AT_FRAC,CUR_NODAL_PRESS_AT_FRAC,STEP,Sigma_XX_INIT,Sigma_YY_INIT,REF_LEL_FRAC_ELEM, &
& FRACTURED_DISP_LOAD,M2_Matr,Prev_P_P,Permiab,Mue,CT,Phi,ALPHA,D,Curr_Pw,Prev_S_Prop_Pw,PORO_DISPL_BNDY_LOAD,&
& Ux_INIT,Uy_INIT,UX_WITH_EXACT_FRAC_PRESS,Uy_WITH_EXACT_FRAC_PRESS)
PRINT*,"AFTER EXACT MODEL"

DO k=1,nN_X
WIDTH_WITH_EXACT_PRESS(k)=Uy_WITH_EXACT_FRAC_PRESS(FRAC_X_AXIS_NODE(k))
END DO

!AFTER THE DISPLACEMENT CONVERGENCE
!refNode_Upd=0.0
refNode_Upd(:,1)=refNode(:,1)-Ux_WITH_EXACT_FRAC_PRESS
refNode_Upd(:,2)=refNode(:,2)-Uy_WITH_EXACT_FRAC_PRESS

! START THE CONVERGENCE OF CRITICAL WIDTH
UY_FREE_NODE_AT_X(1)=Uy_WITH_EXACT_FRAC_PRESS(REF_LEL_FRAC_ELEM(i,2))
UY_FREE_NODE_AT_X(2)=Uy_WITH_EXACT_FRAC_PRESS(REF_LEL_FRAC_ELEM(i,3))


WIDTH_I_FRAC=ABS(UY_FREE_NODE_AT_X(1)*2.0)
PRINT*,"WIDTH_I_FRAC",WIDTH_I_FRAC
DIFF_WIDTH=ABS(WIDTH_CRIC_FRAC-WIDTH_I_FRAC)
TOLLARENCE=DIFF_WIDTH/WIDTH_CRIC_FRAC

!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
IF (TOLLARENCE.GT.TOLLARENCE_ACC) THEN
CALL CRITC_WIDTH_COVG(WIDTH_CRIC_FRAC,WIDTH_I_FRAC,Curr_Pw,MinHStr,G_PRESS,LOOP,FINAL_Pw)
Curr_Pw=FINAL_Pw
FINAL_Pw=0.0
LOOP=LOOP+1
END IF
END DO ITERATIVE ! ITERATION END FOR CRITICAL

 

[mathmate]

It's a little far-fetched, but can you try reversing the order of the allocation of memory, in case the memory was fragmented?

Change from
ALLOCATE(TotSol(nN*3),M_Load(nN*3),BANDED_PORO_MAS S_MAT(2*BAND_WIDTH+1,nN*3))

to:
ALLOCATE(BANDED_PORO_MAS S_MAT(2*BAND_WIDTH+1,nN*3),TotSol(nN*3),M_Load(nN*3))

If it doesn't work, try separate the allocation statement into three and print statements like:
before allocation of matrix
after allocation of matrix
before allocation of vector TotSol
after...

 

[Rahman_M_M]

Thanks mathmate,

I follow your instruction and find that before allocation of vector M_Load i get that message and program terminates. That means program is unable to allocate that vector(M_Load). Anything particular I can do to solve this problem now?

Thanks again
Rahman_M_M

 

[mathmate]

I would try a much smaller problem.

If it gives the same problem, we'll have to look at the programming aspects, such as the correct memory allocation routines, parameters, etc.

If the smaller program works, it means that you are lacking memory. You may want to see if there are configuration parameters that allow you to increase the useable memory. Note that the fact that you have 3Gb of memory on your computer does not automatically mean that it will be all used for your fortran program. Can you calculate how much memory your program requires for the program that does fails? Is this amount commensurate with the memory available?

 

[Rahman_M_M]

Thanks mathmate for your quick response. I try with smaller problem it doesn't work. By the by in our school other student use more bigger problem( They handle 9000x9000 matrix) with same configuration. So I think it's not lacking memory problem.


Regards.
Rahman_M_M

 

[mathmate]

If I remember right, the programme fails on the second run.
Is it still the case with your latest tests, i.e. it works for the first run, but fails on the second?

 

[mathmate]

I must admit that I am not familiar with Fortran 90. According to the reference manual, it does not look like we can allocate many arrays in the same statement, but rather in separate allocate statements just like the way you did with the deallocate statements.
Could you try changing from:
ALLOCATE(TotSol(nN*3),M_Load(nN*3),BANDED_PORO_MAS S_MAT(2*BAND_WIDTH+1,nN*3))
to
ALLOCATE(TotSol(nN*3))
ALLOCATE(M_Load(nN*3))
ALLOCATE(BANDED_PORO_MAS S_MAT(2*BAND_WIDTH+1,nN*3))

 

[Rahman_M_M]

I tried with separate allocate statements. Still it 's not working and it fails in second run.

 

[mathmate]

If the allocation fails on the second attempt, you'd have to make sure the previous allocation of memory had been deallocated. You can confirm this by a few print statements similar to what you did for the allocate statement.

If allocate still fails at the second attempt, you can add parameters
stat=nStat, errmsg=cErrmsg
which will hopefully tell you what went wrong.
nStat=
0 No error
1 Error in system routine attempting to do allocation
2 An invalid data object has been specified for allocation
3 Both error conditions 1 and 2 have occurred
and cErrmsg can be a character(32) variable.

You can put these two parameters in the deallocate statements as well.

 

[mathmate]

Just found something suspicious:

 ALLOCATE(BANDED_PORO_MAS S_MAT(2*BAND_WIDTH+1,nN*3))

has an extra space in the word MASS.
Could you check if that caused the problem, since it may be interpreted as a different variable name?

 

[Rahman_M_M]

In actual code no space between the word mass. I already deallocate the matrix. If we do not deallocate a variable and try to allocate again the program gives us the error message variable is already allocated.

 

[mathmate]

If we do not deallocate a variable and try to allocate again the program gives us the error message variable is already allocated.

This is why I want to double-check, since you mentioned that it is the second run that fails.
For example, do not deallocate and see if it DOES give a message, or use the status parameters to check that the deallocation has actually taken place.

 

[Rahman_M_M]

So it's not a deallocation problem. I try in different way, I declare the variable BANDED_PORO_MAS S_MAT(2*BAND_WIDTH+1,nN*3) in normal way, not allocatable. still i get same problem.

 

[mathmate]

If you declare size parameters in a subroutine call, and use them to declare variables inside of the subroutine, the memory will be allocated on entry, and deallocated on exit.
Try a test of the form:
SUBROUTINE CALC(N1,N2)
REAL A(N1),B(N2)
...
RETURN
END
See if it works. Then you do not have to worry about allocation or deallocation.

 

[mathmate]

Do you get the 'same' problem meaning that it fails on the second run?

Back to the original code: if you put the stat=nStat parameter in the allocate and (especially) deallocate calls, what codes did you get in the first run, the second?

I assume you need to use the allocate statement because the size varies from run to run, since the bandwidth usually varies with different problems. If you have time, you can try the object time dimensioning I suggest in my previous post using a simple subroutine. Unfortunately I do not have access to F90 to test it before. If it works, then you can try incorporating the idea in your program. It is supposed to automatically allocate and deallocate, but the size parameters must be formal parameters.