problem with 2d arrays in fortran 90!

problem with 2d arrays in fortran 90!

I have written a simple program,which should transform 1d array into 2d.But I get something stupid.I know that FORTRAN is writing arrays in columns,so just to ask is it compiler issue.I am using Intel Fortan Composerxe 2011.
The program:
program pr2
implicit none
integer:: i,j
real, dimension(21,80) :: istat
open(7,file='sens.dat',status='old',action='read')
do i=1,21
read(7,*)(istat(i,j),j=1,80)
write(*,*)istat
end do
end

I am expecting 21times 80 table.But:
9.0841827E-16 -1.1127702E-18 -3.7055020E-18 -1.7058567E-17 -2.4920699E-17
-3.1891962E-17 -1.1954544E-21 1.2890473E-14 -1.0779296E-21 9.0450016E-18
6.0385893E-18 -1.2223998E-16 3.2401010E-20 -6.3506461E-18 -5.0402904E-17
-6.1655475E-10 4.2683850E-23 -1.0390595E-16 1.2178668E-15 -1.7920678E-12
1.2697568E-08 -3.8678227E-05 1.7134852E-08 -2.8258935E-08 -4.4748720E-08
-3.2356020E-08 -2.4845672E-08 -2.2212359E-11 -2.2761658E-04 -6.7894028E-12
-9.0211152E-08 -1.6370835E-07 -1.7512080E-07 -7.6558049E-10 -4.1889128E-09
2.2164288E-08 -2.586304 8.4790370E-14 -2.4060276E-07 3.0392204E-05
I have attached input file.

AttachmentSize
Downloadchemical/x-mopac-input sens.dat26.25 KB
19 posts / 0 new
Last post
For more complete information about compiler optimizations, see our Optimization Notice.

No, it is not a compiler issue, nor is the result that you received unexpected.

An array in an I/O list is printed out in the order Col-1, Col-2, ... Secondly, when you do not specify a format, the compiler prints the output using a default format that may vary from one compiler to another. Thirdly, since your WRITE statement is inside the loop, the whole array is printed out 21 times; the first 20 times, the contents of the array have not been fully initialized, so many of the values printed out are whatever junk happened to be in memory.

Your write statement is within the do...enddo construct, which means that the whole istat array will be written 21 times. Placing the write statement after the enddo should solve the problem. One alternative (probably the one you want) is to print just the array slice you just read, as in:

program pr2
    implicit none
    integer:: i,j
    real, dimension(21,80) :: istat
    open(7,file='sens.dat',status='old',action='read')
        do i=1,21
            read(7,*) istat(i,:)
            write(*,*) istat(i,:)
            read *  !this is just a PAUSE-like statement
        end do
    close (7)
end program pr2

You might also want to explicitly state the desired format, instead of using *.

Well thanks guys but even with John,s code I do not get what I need.I NEED 21x80 array which has physical meaning,then I would take transpose of it and so on...
So i need somehow to specify fomat to write me this.

You are missing the obvious: you have just 1 value per line. So your implied DO above for j will not work. You need this:

do j=1,80
do i=1,21
read(7,*)istat(i,j)
end do
end do

Read 1 item per line, read 21 elements down a column along the row index, do this for each of the 80 columns.

ron

Ron,it is still working the same.
program pr2
implicit none
integer:: i,j
real, dimension(21,80) :: istat
open(7,file='sens.dat',status='old',action='read')
do j=1,80
do i=1,21
read(7,*)istat(i,j)
end do
end do
write(*,*)istat
end
I got five column output.
4.8205439E-13 1.4917772E-03 -3.5040978E-02 2.0828361E-04 -5.2204872E-12
1.9907949E-02 -3.0700287E-02 1.6273017E-04 -1.8292824E-12 6.8864161E-03
-2.1972844E-02 6.4284824E-05 -5.2571250E-13 1.9359445E-03 -1.3504234E-02
3.6830806E-05 -9.6655982E-12 3.5901170E-02 2.5154045E-04 -1.5395350E-05
-6.2947226E-14 3.5201247E-05 -1.8068762E-02 -1.1046973E-04 -6.2351074E-13

If you want as specific output format, specify it yourself. list-directed format is chosen by each library, without knowledge of what is in your source code or data.

exactly, what is it you are looking for in the output? Do you want 21 lines with 80 values in each line? Do you want it comma separated so you can pull it into a spreadsheet or something. You are not being clear in your expectations.

ron

Ok yes It works fine with proper formating.

Ron, he is using list-directed input, so record terminators in the input are equivalent to blanks separating the numbers, and John's modified version will work fine -- as long as the input file has 21 X 80 numbers and does not contain slashes.

I have problems again.
program pr

c creat 2d array appres and phase
implicit none
integer::a,b
real,dimension(20,4)::model,data
open(20,file='apandph.dat',status='unknown')
do a=1,20
read(20,175)(model(a,b),b=1,4)
175 format(f10.6,2x,f11.8,2x,f10.6,2x,f11.8)
write(*,*)(model(a,b),b=1,4)
end do
end
This works fine.But in the main program it can not open.My strategy is to write it as 2d array,later to manipulate them.
forrtl: severe (24): end-of-file during read, unit 20, file /home/milenko/ircg/apandph.dat
Image PC Routine Line Source
ircg 0809B103 Unknown Unknown Unknown
ircg 08099E20 Unknown Unknown Unknown
Relevant part of the code:

ro_e=k*(zxy_re_e**2+zxy_im_e**2)
phase_e=atan(zxy_im_e/zxy_re_e)
ro_h=k*(zxy_re_h**2+zxy_im_h**2)
phase_h=atan(zxy_im_h/zxy_re_h)
write(20,175)ro_e,phase_e,ro_h,phase_h
175 format(f10.6,2x,f11.8,2x,f10.6,2x,f11.8)
end do

c creat 2d array appres and phase
c do b=1,4
do a=1,20
read(20,*)(model(a,b),b=1,4)
c175 format(f10.6,2x,f11.8,2x,f10.6,2x,f11.8)
end do
Everything is the same and it is NOT WORKING.

When I try to read data:
milenko@milenkons:~/ircg1$ ./ircg
20
-0.4799119 3.563040 -0.5209458 0.0000000E+00

data.dat
3.811112 -0.47991191 3.563044 -0.52094584
3.489387 -0.48291942 3.591924 -0.53019292
3.600001 -0.47000002 3.742345 -0.58005581
9.501012 -0.33234569 6.500099 -0.66299910
125.008721 -0.97001245 202.231452 -0.77712679
108.002329 -0.98238911 191.123421 -0.99128101
54.856758 -0.72918271 53.887214 -1.01232191
36.189101 -0.46877231 8.612981 -0.83198211
88.000000 -0.81200501 585.000000 -0.62143009
79.000000 -0.80809071 81.200000 -0.70503449
57.524508 -0.72580025 9.512349 -0.89342348
53.210210 -0.81432536 11.000000 -0.72104921
11.545436 -0.39213345 3.9876001 -1.21399992
72.000000 -1.03245234 204.008761 -1.10023459
100.503948 -1.22038476 154.238451 -1.12498374
33.293858 -1.15092732 37.000873 -1.19872445
43.988878 -1.24325900 55.905234 -1.16998874
108.954144 -0.68512361 103.232623 -0.67260660
108.114144 -0.65512361 103.232623 -0.64260660
108.018521 -0.61990360 110.122787 -0.62221227

Doesn't see the first one.

In your last post, you did not show the declarations of the variables in the I/O list, the READ statement and the FORMAT used. Nor did you provide the WRITE statements that gave the first two lines of output that you showed.

With such scanty information, I cannot help.

Ok,here is the code:
module constants
c ================
c Often used constants in the EM modeling code
implicit none
c
c Machine precision for REAL*4 values
real*4 mach_real4
parameter(mach_real4=tiny(mach_real4))
c
c Imaginary unit and PI
complex*8 imc
real*4 pi
parameter(pi=3.141592653589793,imc=(0.,1.))
c
c Scale factors used in re-calculating angle values from degrees
c into radians and vice versa
real*4 prev_rad2deg,prev_deg2rad
parameter(prev_rad2deg=180./pi,prev_deg2rad=pi/180.)
c
c Scale factors used in changing the units of the impedance; standardly
c SI units (ohm) are used; if practical units (mV/km/nT) are preferred
c then uncomment the second parameter PREV_ZUNIT definition
real*4 prev_zunit,prev_siunit,prev_pracunit
parameter(prev_siunit=1.,prev_pracunit=2500./pi)
c
c Natural logarithm of 10
real*4 aln10
parameter(aln10=alog(10.))
c
c Vacuum permeability (used throughout the code)
real*4 mu0
parameter(mu0=4.e-7*pi)
end module

module params
c =============
c Maximum values of parameters of the 2D model
implicit none
c
c Maximum number of mesh lines (including the margins) in the
c horizontal and vertical direction
integer*4 nmax,mmax
parameter(nmax=41,mmax=35)
c
c Maximum number of periods and of various MT transfer functions
integer*4 npmax,nfmax
parameter(npmax=20,nfmax=4)
c
c Maximum number of uniform conductivity domains in the 2D model
integer*4 ncmax
parameter(ncmax=(nmax-1)*(mmax-1))
end module

module mt2dmod
c ==============
c Parameters of the 2D model
use params
c
implicit none
c
c Horizontal and vertical mesh steps and horizontal and vertical
c coordinates of the mesh lines (including the margins of the model)
real*4 sy(nmax-1),sz(mmax-1),syy(nmax),szz(mmax)
c

integer*4 n,m,ma,nc,nbl,nbr
c
c Values of resistivities and conductivities in the 2D model
real*4 res(ncmax),cond(ncmax)
c
c Flags indicating that a sensitivity with respect to the particular
c conductivity is to be computed (1) or not (0)
integer*4 ivar(ncmax)
c
c Conductivity map of the 2D model; each cell is attributed the
c conductivity COND according to the index IC
integer*4 ic(mmax-1,nmax-1)
c
c Thicknesses and conductivities of the 1D marginal models at
c "infinity" for the left and right hand side margin of the 2D model
real*4 hbl(mmax),sbl(mmax),hbr(mmax),sbr(mmax)
c
c Conductivity indices of the layers of the 1D marginal models given
c in terms of the indices of the conductivity map IC
integer*4 isbl(mmax),isbr(mmax)
c
c Boundary conditions for the 2D model at the margings, Left, Right,
c Top, and Bottom. The first index is for the polarization (1 for E,
c 2 for H)
complex*8 bol(2,mmax),bor(2,mmax),bot(2,nmax),bob(2,nmax)
end module
c
program inversion regularized conjugate gradient

use constants
use params
use mt2dmod

implicit none

integer*4 n1,i,j,icd,ivarprd
integer*4 iprd,iprd0,isite,ico
integer*4 np,nsites
real*4 per(npmax)
integer*4 isit(nmax),ircpr,m1
integer*4 ivarp(ncmax)
integer:: a,b
real,dimension(20,4)::model,data,dmf
real:: period,d,d1,g,g1,k
real:: zxy_re_e,zxy_im_e,zxy_re_h,zxy_im_h
real:: sens_zxy_ree,sens_zxy_ime,sens_zxy_reh,sens_zxy_imh
real:: ro_e,ro_h,phase_e,phase_h
real:: sens_ro_e,sens_phase_e,sens_ro_h,sens_phase_h
real:: ro_ed, phase_ed,ro_hd,phase_hd
real:: dmro_e,dmphase_e,dmro_h,dmphase_h

c open periods and isit(i)s of interest
open(2,file='test_p7.dat',status='old')
c read conductivity map
open(3,file='test_sh7.dat',status='old')
c read direct and sensitivity from mt2d(results from josef code)
open(10,file='test_ie.dat',status='old')
open(11,file='test_sense.dat',status='old')
open(12,file='test_ih.dat',status='old')
open(13,file='test_sensh.dat',status='old')
open(14,file='sens_results.dat',status='unknown')
open(15,file='data.dat',status='old')
c open(16,file='datamisfit.dat',status='unknown')
c open(17,file='sens_transpose.dat',status='unknown')
open(20,file='apandph.dat',status='unknown')

read(3,*)n,m,ma
n1=n-1
m1=m-1
c
c Read the horizontal and vertical mesh steps, in km
read(3,*)(sy(i),i=1,n1)
read(3,*)(sz(i),i=1,m1)
syy(1)=0.
do i=1,n1
syy(i+1)=syy(i)+sy(i)
enddo
szz(ma)=0.
do i=ma,m1
szz(i+1)=szz(i)+sz(i)
enddo
do i=ma-1,1,-1
szz(i)=szz(i+1)-sz(i)
enddo
do i=1,ma-1
do j=1,n1
ic(i,j)=1
enddo
enddo

read(3,*)(ic(i,j),j=1,n1)
enddo
icd=1
res(icd)=-1.
cond(icd)=0.
ivar(icd)=0
c
c Read the actual number of different conductivities used in the conductivity
c map. This number is for conductivities in the Earth only, the air domain is
c NOT included here!
read(3,*)nc
nc=nc+1
ivarprd=0

do i=2,nc
read(3,*)icd,res(icd),ivar(icd)
cond(icd)=1./res(icd)
if(ivar(icd).ne.0)then
ivarprd=ivarprd+1
ivarp(ivarprd)=icd
endif
enddo
c
c That is all to be read from the structure file. Check especially the
c correspondence between the indices in the conductivity map and those
c in the resistivity list. This check is NOT carried out in the code!
close(3)
c
c
c Inputs from the period/isit(i)s file.
c Actual number of periods is read and a list of periods in seconds
read(2,*)np
read(2,*)(per(i),i=1,np)

read(2,*)nsites
read(2,*)(isit(i),i=1,nsites)
write(*,*)nsites
c
c Flag whether the reciprocity approach is to be used for the sensitivity
c calculations (IRCPR=1) or not (IRCPR=0)
read(2,*)ircpr
c
c That is all the input from the period/structure file
close(2)
c Read reslults from Josef's code
do n=1,nsites
isite=isit(i)
read(10,150)period,isit(i),iprd0,zxy_re_e,zxy_im_e
150 format(f12.4,2i5,2e15.6)
c write(*,*)period,isit(i),iprd0,zxy_re_e,zxy_im_e
read(12,152)period,isit(i),iprd0,zxy_re_h,zxy_im_h
152 format(f12.4,2i5,2e15.6)
c write(*,*)period,isit(i),iprd0,zxy_re_h,zxy_im_h
k=period/(2*mu0)
c write(*,*)k
ro_e=k*(zxy_re_e**2+zxy_im_e**2)
phase_e=atan(zxy_im_e/zxy_re_e)
ro_h=k*(zxy_re_h**2+zxy_im_h**2)
phase_h=atan(zxy_im_h/zxy_re_h)
write(20,175)ro_e,phase_e,ro_h,phase_h
175 format(f10.6,2x,f11.8,2x,f10.6,2x,f11.8)
c write(*,*)ro_e,phase_e,ro_h,phase_h
end do

c creat 2d array appres and phase
do a=1,20
read(20,175)(model(a,b),b=1,4)
end do
write(*,*)(model(a,b),b=1,4)
c create 2d array data
do a=1,20
read(15,175)(data(a,b),b=1,4)
end do
write(*,*)(data(a,b),b=1,4)

c data misfit.deduct one array from the other
c dmf(a,b)=model(a,b)-data(a,b)
c write(16,*)dmf(a,b)

c sensitivity read and calculations

do n=1,nsites
isite=isit(i)
do ico=1,nc
if(ivarp(ico).ne.0)then
iprd=ivarp(ico)
read(11,151)period,isit(i),iprd,sens_zxy_ree,sens_zxy_ime
151 format(f12.4,2i5,2e15.6)
c write(*,*)period,isit(i),iprd,sens_zxy_ree,sens_zxy_ime
read(13,153)period,isit(i),iprd,sens_zxy_reh,sens_zxy_imh
153 format(f12.4,2i5,2e15.6)
c write(*,*)period,isit(i),iprd,sens_zxy_reh,sens_zxy_imh
d=zxy_re_e**2+zxy_im_e**2
d1=1/d
g=zxy_re_h**2+zxy_im_h**2
g1=1/g
sens_ro_e=2*k*(sens_zxy_ree*zxy_re_e+sens_zxy_ime*zxy_im_e)
sens_phase_e=d1*(sens_zxy_ime*zxy_re_e-zxy_im_e*sens_zxy_ree)
sens_ro_h=2*g*(sens_zxy_reh*zxy_re_h+sens_zxy_imh*zxy_im_h)
sens_phase_h=g1*(sens_zxy_imh*zxy_re_h-zxy_im_h*sens_zxy_reh)
write(14,*)sens_ro_e
write(14,*)sens_phase_e
write(14,*)sens_ro_h
write(14,*)sens_phase_h

end if
end do
end do

end

I am attaching the input files.

Attachments: 

This code contains numerous syntax errors, whether compiled as a .f file or a .f90 file. Line 158 contains what looks like part of a DO loop but there is no matching DO statement preceding it.

You could not have reached the point of seeing run-time errors with code that the compiler rejects. It is a waste of people's time when you ask for help after posting code that has syntax errors and a description of run-time errors that the posted source code could not have produced.

Please post compilable source code using the syntax highlighter (The yellow pencil icon), in a suitable Fortran mode.

Thanks a lot mecej4.I have more versions probably I have put the one that can not be compiled.I am attaching one that is compiled but can not read files 15 and 20.

Attachments: 

AttachmentSize
Downloadapplication/octet-stream ircg.for10.44 KB

There are two main errors in your last program (ircg.for, attached to #16) that are related to file processing, near lines 249-260:

(i) After writing data to unit-20, you proceed to read from the same file. If you intend to read the data that was just written, you must issue a REWIND(20) statement before reading, and it would be useful to check for EOF, ERR in the subsequent READ statement.

(ii) Lines 251-254 as you have them:

[fxfortran]         do a=1,20
read(20,175)(model(a,b),b=1,4)
end do
write(*,*)(model(a,b),b=1,4)
[/fxfortran]

have the problem that the WRITE statement should have been inside the DO loop. The loop counter 'a' has the value 21 after the loop is finished, and model(21,1:4) contains undefined values. Similarly for Lines 256-260.

I observe from your posts that you are not practicing defensive programming. Until you reach competence in the programming language you should compile code that you write with the "-traceback", "-check all" and "-warn all" options. Try compiling your program with these options, before you fix the errors that I pointed out, and run. Learn to use the runtime error messages to fix your source code and see if you can identify the same errors as I did.

Thanks a lot,works fine.

Leave a Comment

Please sign in to add a comment. Not a member? Join today