reading complex data using implied do loops
Ragu
program cread
use iso_fortran_env, only: output_unit
implicit none
integer, parameter :: sp_k = kind(1.0) ! Default Type of Real
integer, parameter :: dp_k = selected_real_kind(2*precision
(1.0_sp_k))
complex (kind = dp_k), dimension(4) :: cdata
real(kind = dp_k), dimension(4) :: redata, imdata
integer :: ii
continue
! Part 1. Code for reading both real and imaginary in an implied loop
open(unit = 11, file = 'data.inp', form = 'formatted', status = 'old')
read(11,'(2E16.6,16X)') (cdata(ii), ii = 1,4)
close(11)
write(output_unit,'(3A16)') 'REAL PART', 'IMAG PART', 'AMPLITUDE'
write(output_unit,'(2E16.6)') (cdata(ii), ii = 1, 4)
! Part 2. Code to read only the imaginary part of complex
open(unit = 11, file = 'data.inp', form = 'formatted', status = 'old')
read(11,'(16X,E16.6,16X)') (imdata(ii), ii = 1,4)
close(11)
redata = 0.0_dp_k
cdata = cmplx(0.0_dp_k,imdata,kind = dp_k)
write(output_unit,'(3A16)') 'REAL PART', 'IMAG PART', 'AMPLITUDE'
write(output_unit,'(2E16.6)') (cdata(ii), ii = 1, 4)
end program cread
data file 'data.inp'
0.100034E+01 -0.382548E-04 0.100034E+01
0.100054E+01 -0.671324E-04 0.100054E+01
0.100077E+01 -0.102652E-03 0.100077E+01
0.100106E+01 -0.144950E-03 0.100106E+01
I am trying to upgrade my code to accomodate complex values.
Previously I dealt only with reals and this is my first day using
complex. My current code is structured and uses a lot of implied do
loops to read data from both formatted and unformatted data. I am
playing with a test code. I am trying to understand some basic methods
for reading complex values using implied do loops. I am not getting a
much help from my fortran references.
I have two questions:
1. How do I print the ABS(cdata) in the same implied do loop ? Or
would you advise me to just print real(cdata(ii)), imag(cdata(ii)), abs
(cdata(ii)) in a normal do loop?
2. Is it possible to avoid the use of a temporary variable in Part 2
using implied do loop ? My actual arrays are huge in size.
I am hoping that the answers to the above will help me address most of
the special cases that might arise. Thanks.
robin
Just put ABS(cdata) in the same implied-DO.
| Or
| would you advise me to just print real(cdata(ii)), imag(cdata(ii)), abs
| (cdata(ii)) in a normal do loop?
| 2. Is it possible to avoid the use of a temporary variable in Part 2
| using implied do loop ?
Yes.
Ragu
> | 1. How do I print the ABS(cdata) in the same implied do loop ?
>
> Just put ABS(cdata) in the same implied-DO.
If it is not a big thing can you please give me a sample syntax. The
combination that I tried showed up as errors.
>
> | Or
> | would you advise me to just print real(cdata(ii)), imag(cdata(ii)), abs
> | (cdata(ii)) in a normal do loop?
> | 2. Is it possible to avoid the use of a temporary variable in Part 2
> | using implied do loop ?
>
> Yes.
>
Umm. Can you provide some details ? Does it involves transfer and
reshaping that I have totally avoided till date ?
m_b_metcalf
Does this solve both your problems?
Regards,
Mike Metcalf
program cread
use iso_fortran_env, only: output_unit
implicit none
integer, parameter :: sp_k = kind(1.0) ! Default Type of Real
integer, parameter :: dp_k = selected_real_kind(2*precision
(1.0_sp_k))
complex (kind = dp_k), dimension(4) :: cdata
real (kind = dp_k), dimension(2*4) :: tdata = 0.0_dp_k
equivalence (tdata, cdata)
integer :: ii
! Part 1. Code for reading both real and imaginary in an implied loop
open(unit = 11, file = 'data.txt', form = 'formatted', status =
'old')
read(11,'(2E16.6,16X)') (cdata(ii), ii = 1,4)
write(output_unit,'(3A16)') 'REAL PART', 'IMAG PART', 'AMPLITUDE'
write(output_unit,'(3E16.6)') (cdata(ii),abs(cdata(ii)), ii = 1, 4)
! Part 2. Code to read only the imaginary part of complex
rewind 11
cdata = 0
read(11,'(16X,E16.6,16X)') (tdata(ii), ii = 2,2*4, 2)
close(11)
write(output_unit,'(3A16)') 'REAL PART', 'IMAG PART', 'AMPLITUDE'
write(output_unit,'(3E16.6)') (cdata(ii),abs(cdata(ii)), ii = 1, 4)
end program cread
Ragu
careful with the backets and look dumb after shown the way.
The second part deals uses equivalence. Something I have read a lot in
this forum but never used it. So having the equivalence of two arrays
of same data type, kind and length shares the same storage unit and
same memory location ? The book talks about storage unit but not
memory.
The arrays in my real code are all allocatable arrays. Don't I need to
allocate for the tdata ? In that case will I not be using extra
memory? Is there any way to avoid such situation?
Thanks. I guess many stayed out of this thread thinking that it is
home work.
m_b_metcalf
You can't use equivalence with allocatable arrays (MR&C, Section
20.2.2), so you'll have to cheat another way. The subroutine has to be
externaal. It still may require switching off interface checks!
Regards,
Mike Metcalf
! Part 2. Code to read only the imaginary part of complex
rewind 11
cdata = 0
call read_im(cdata)
close(11)
write(output_unit,'(3A16)') 'REAL PART', 'IMAG PART', 'AMPLITUDE'
write(output_unit,'(3E16.6)') (cdata(ii),abs(cdata(ii)), ii = 1, 4)
end program cread
subroutine read_im(im)
integer, parameter :: sp_k = kind(1.0) ! Default Type of Real
integer, parameter :: dp_k = selected_real_kind(2*precision
(1.0_sp_k))
real(kind = dp_k), dimension(2, 4) :: im
read(11,'(16X,E16.6,16X)') (im(2, ii ), ii = 1, 4)
end subroutine read_im
Richard Maine
> The second part deals uses equivalence. Something I have read a lot in
> this forum but never used it. So having the equivalence of two arrays
> of same data type, kind and length shares the same storage unit and
> same memory location ? The book talks about storage unit but not
> memory.
For practical purposes, storage and memory are just different words for
the same thing. There is a subtle theoretical distinction in that
"storage unit" is a term used by the standard, while "memory" could be
called an implementation. However, for your purposes, this distinction
is unimportant.
> The arrays in my real code are all allocatable arrays.
Equivalence doesn't work with allocatable.
I have trouble comming up with a method for reading just the imaginary
part of a complex array that is both simple and standard conforming.
There are certainly hacks that would probably work, but I'll not
recommend them. There have been proposals to add a component-like syntax
for the real and imaginary parts of a complex. I support the idea. That
would make things like this simple. However, it isn't in the current
standard. (I don't recall whether that one is in the f2008 draft).
Since your arrays are allocatable, the simplest solution that occurs to
me is to go ahead and use a temporary array. You mention that the arrays
are "huge". That is some argument against it, but I'm not sure it is an
overwhelming one.
If you just allocate and deallocate the temporary array right around the
READ statement, the memory it uses will be available for other purposes
at other times in the program execution. The I/O parts of most
scientific programs are usually not the most memory intensive, so there
is a good chance that using the temporary array there won't significan't
increase your overall memory usage. It could happen, but I wouldn't
assume so without more information.
The time to do the copy from the temporary array is also unlikely to be
large relative to other things. If nothing else, I/O is slow. Just the
time to do the I/O would probably dwarf the time of a single copy,
without even counting anything else that you might later do with the
array.
I'm not an expert in I/O implementation, but it would not surprise me if
the I/O system didn't make a temporary array in its implementation of
reading into a non-contiguous array anyway. (The imaginary part of a
complex array would be noncontiguous). So going to a lot of trouble to
take such a temporary array out of the visible code, might just have the
effect of moving it behind your back.
--
Richard Maine | Good judgment comes from experience;
email: last name at domain . net | experience comes from bad judgment.
domain: summertriangle | -- Mark Twain
Richard Maine
> You can't use equivalence with allocatable arrays (MR&C, Section
> 20.2.2), so you'll have to cheat another way. The subroutine has to be
> externaal. It still may require switching off interface checks!
That's exactly the hack that was foremost in my mind when I said (in the
post that I was apparently composing at the sam etime that Mike was
doing his)
>> There are certainly hacks that would probably work, but I'll not
>> recommend them.
--
James Van Buskirk
news:1jd0vbn.4ahxp51iptv22N%nos...@see.signature...
> I have trouble comming up with a method for reading just the imaginary
> part of a complex array that is both simple and standard conforming.
> There are certainly hacks that would probably work, but I'll not
> recommend them. There have been proposals to add a component-like syntax
> for the real and imaginary parts of a complex. I support the idea. That
> would make things like this simple. However, it isn't in the current
> standard. (I don't recall whether that one is in the f2008 draft).
f2008 has z%RE and z%IM but these aren't implemented yet in gfortran.
You can hack with C_PTR, though:
C:\gfortran\clf\read_imag>type data.inp
0.100034E+01 -0.382548E-04 0.100034E+01
0.100054E+01 -0.671324E-04 0.100054E+01
0.100077E+01 -0.102652E-03 0.100077E+01
0.100106E+01 -0.144950E-03 0.100106E+01
C:\gfortran\clf\read_imag>type read_imag.f90
program read_imag
use ISO_C_BINDING, only: C_PTR, C_LOC, C_F_POINTER
use ISO_FORTRAN_ENV, only: output_unit
implicit none
integer, parameter :: dp = kind(1.0d0)
complex(dp), allocatable, target :: cdata(:)
type(C_PTR) pc
real(dp), pointer :: rdata(:)
integer ii
allocate(cdata(4))
pc = C_LOC(cdata(1))
call C_F_POINTER(pc,rdata,[2*size(cdata)])
open(11,file='data.inp',status='old')
cdata = 0
read(11,'(16X,E16.6,16X)') (rdata(2*ii),ii=1,size(cdata))
close(11)
write(output_unit,'(3A16)') 'REAL PART', 'IMAG PART', 'AMPLITUDE'
write(output_unit,'(3E16.6)') (cdata(ii),abs(cdata(ii)),ii=1,size(cdata))
end program read_imag
C:\gfortran\clf\read_imag>gfortran read_imag.f90 -oread_imag
C:\gfortran\clf\read_imag>read_imag
REAL PART IMAG PART AMPLITUDE
0.000000E+00 -0.382547E-04 0.382547E-04
0.000000E+00 -0.671323E-04 0.671323E-04
0.000000E+00 -0.102652E-03 0.102652E-03
0.000000E+00 -0.144949E-03 0.144949E-03
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
Richard Maine
> "Richard Maine" <nos...@see.signature> wrote in message
> news:1jd0vbn.4ahxp51iptv22N%nos...@see.signature...
>
> > I have trouble comming up with a method for reading just the imaginary
> > part of a complex array that is both simple and standard conforming.
> > There are certainly hacks that would probably work, but I'll not
> > recommend them. There have been proposals to add a component-like syntax
> > for the real and imaginary parts of a complex. I support the idea. That
> > would make things like this simple. However, it isn't in the current
> > standard. (I don't recall whether that one is in the f2008 draft).
>
> f2008 has z%RE and z%IM but these aren't implemented yet in gfortran.
> You can hack with C_PTR, though:
[sample code elided]
Good point.
I should have thought of that one, as I have noted before that C_PTR can
be used for several things unrelated to C, type cheating being on the
list. Unlike type-cheating by taking advantage of argument sequence
association, this one should not need anything to "hide" the trick from
the compiler. It is at least arguably standard conforming. (The
potential arguments are nitpicky enough that I don't feel it worth
spending the time to study, and I'm not even sure which side I'd end up
on; I'll just go with "it ought to work.")
It is still a bit on the hacky side - the kind of code where the intent
is not obvious without comments. I'd probably still argue for making a
temporary array as being more readable, but this might be my second
choice (if I assumed I could limit myself to compilers that had
implemented that part of f2003). My first choice in some future would be
the z%IM thing, but that isn't available yet.
robin
On Jan 27, 6:58 pm, "robin" <robi...@bigpond.com> wrote:
>> | 1. How do I print the ABS(cdata) in the same implied do loop ?
>
>> Just put ABS(cdata) in the same implied-DO.
>If it is not a big thing can you please give me a sample syntax. The
>combination that I tried showed up as errors.
Your implied-DO was: (cdata(ii), ii = 1, 4)
To include ABS(cdata) it becomes:
(cdata(ii), ABS(cdata(ii)), ii = 1, 4)
>> | Or
>> | would you advise me to just print real(cdata(ii)), imag(cdata(ii)), abs
>> | (cdata(ii)) in a normal do loop?
>> | 2. Is it possible to avoid the use of a temporary variable in Part 2
>> | using implied do loop ?
>
>> Yes.
>Umm. Can you provide some details ? Does it involves transfer and
>reshaping that I have totally avoided till date ?
No.