forums.silverfrost.com Forum Index forums.silverfrost.com
Welcome to the Silverfrost forums
 
 FAQFAQ   SearchSearch   MemberlistMemberlist   UsergroupsUsergroups   RegisterRegister 
 ProfileProfile   Log in to check your private messagesLog in to check your private messages   Log inLog in 

error in array construction

 
Post new topic   Reply to topic    forums.silverfrost.com Forum Index -> Support
View previous topic :: View next topic  
Author Message
AntonioGJerez



Joined: 01 Jun 2015
Posts: 9
Location: Spain

PostPosted: Mon Jun 01, 2015 1:51 pm    Post subject: error in array construction Reply with quote

Dear Friends

This code produces a runtime Access Violation Error when trying to construct an array which includes elements from an array-valued function

Moreover, If I define REAL kk(5), then kk=(/ 0.0, diag4(matrix) /) gives a compilation error "Non-conformant array shapes infirst Rank of an array expression (5 and 2)"
It seems that diag4(matrix) is interpreted as a real instead of as a 4 elements array.

Any sugestion?

Thank you in advance

Antonio


program bug

interface diag4
! Extracts the diagonal from a 4x4 matrix (assumed to be square)
function diag4(A)
real,dimension(4,4),intent(in)::A
real,dimension(4)::diag4
end function diag4
end interface

real matrix(4,4)
matrix=reshape((/1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,11.,12.,13.,14.,15.,16./),(/4,4/))
print*,diag4(matrix) ! this runs ok
print*,(/ 0.0, diag4(matrix) /) ! this leads to Access Violation Error
end program

function diag4(A)
implicit none
real,dimension(4,4),intent(in)::A
real,dimension(4)::diag4
integer index
do index=1,size(diag4)
diag4(index)=A(index,index)
enddo
return
end function
Back to top
View user's profile Send private message
JohnCampbell



Joined: 16 Feb 2006
Posts: 2239
Location: Sydney

PostPosted: Tue Jun 02, 2015 8:34 am    Post subject: Reply with quote

I find the following works:
real out(5)
...
out = 0
out(2:5) = diag4(matrix)
print*, out

I am not sure if out = (/ 0., diag4(matrix) /) should work, but I suspect not, although
out = (/ 0., (/ 1., 6., 11., 16./) /) or
out = (/ (/ 0., 1./), (/ 6., 11., 16./) /)
does appear to work.

John
Back to top
View user's profile Send private message
AntonioGJerez



Joined: 01 Jun 2015
Posts: 9
Location: Spain

PostPosted: Tue Jun 02, 2015 8:00 pm    Post subject: Reply with quote

Hi John

thanks, nevertheless, I would like to include the array (/ 0., diag4(matrix) /) into a more complicated expression,

e.g. ... * (/ 0., diag4(matrix) /)+...

it should be better to avoid storing (/ 0., diag4(matrix) /) in a variable in several steps.

I have checked that the above array constructor does work a scalar function defined as

function diag1(A)
real diag1
...

nevertheless it does not work indeed with a one element array:

function diag1(A)
real,dimension(1)::diag1
...

regards
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1429

PostPosted: Wed Jun 03, 2015 11:31 pm    Post subject: Reply with quote

Your post touches multiple issues, and I think that we can tackle them better if we take them one at a time.

1. Matlab and Fortran are different. Matlab is an interpreted language and specialized for one thing: manipulating matrices. It has rather weak notions of types and is not suited for general computing. Fortran is a compiled language and newer versions of Fortran (e.g. F95) have stricter rules regarding types and mixing types. A discussion of the Matlab way vis-a-vis the Fortran way is probably not appropriate in the Silverfrost forums.

2. When the compiler rejects a piece of code, we may have to consult the relevant Fortran standard to interpret the error messages given by the compiler. If the reading of the standard indicates that the code is conforming, file a compiler bug report with the shortest example code that demonstrates the bug, and cite the relevant sections of the standard. Your code demonstrates such a bug, and you may proceed to file a bug report. When you do, refrain from describing alternative codes that work, but do state if you have tried other compilers and they compile the code correctly even though FTN95 does not.

3. If there is a certain idiom that best expresses a calculation, but the language standard in effect (F95 for FTN95) does not allow it, describe the idiom and ask for alternative ways of expressing the idiom in standard Fortran.

The following program falls into the third category that I just described. It uses the "implied DO loop" construct to assemble a 1-D array using one or more constants and the diagonal elements of a 2-D array. Note that it avoids using a function with a return value that is an array. Such functions are problematic because they require an explicit interface wherever they are used and, if no such interface is provided, cause erratic run time behaviour.
Code:
program nobug
implicit none
integer, parameter :: M = 4
real, parameter :: Pi = 3.1415927
real matrix(M,M)
real :: diag(M+1)
integer :: i
matrix=reshape((/ 1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,11.,12.,13.,14.,15.,16. /), &
               (/ M,M /))
diag=(/ Pi, (matrix(i,i), i=1,M) /)
print *,diag
print *,(/ Pi, (matrix(i,i), i=1,M) /)
end program
Back to top
View user's profile Send private message
AntonioGJerez



Joined: 01 Jun 2015
Posts: 9
Location: Spain

PostPosted: Thu Jun 04, 2015 12:48 am    Post subject: Thanks. Nevertheless... Reply with quote

Ok, thanks for your code (even though the proposal of supporting an arbitrary function giving a vector would be much more powerful).

Nevertheless, you say "Such functions are problematic because they require an explicit interface wherever they are used and, if no such interface is provided, cause erratic run time behaviour"

Note that in the original example of this post, the interface is provided, the function works but the array constructor desnīt work

Antonio
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1429

PostPosted: Thu Jun 04, 2015 4:08 am    Post subject: Reply with quote

It is possible to use an implied DO loop to write a function that gives the diagonal of an arbitrary matrix. Here is the code. Note that I had to provide an alternative syntax in Line-6 to work around the FTN95 bug with array constructors.

This code is open to criticism for allocating an array for the function result every time that the function is used. The rules of Fortran 2003 may help to overcome this type of "memory leak" with the "allocate on assignment" feature.

Code:
program getdiag
implicit none
real matrix(2,2),d3(3)
!
matrix=reshape((/ 1.,2.,3.,4. /),(/2,2/))
d3(1) = 0; d3(2:3)=diagf(matrix)
! d3 = (/ 0.0, diagf(matrix) /) !should work, but FTN95 bug is in the way
print*,d3

contains
   function diagf(A)
      real, dimension(:,:), intent(in) :: A
      real, dimension(:),pointer :: diagf
      integer :: i,m
      m=min(size(A,1),size(A,2))
   !
      allocate(diagf(m))
      diagf = (/ (A(i,i), i=1,m) /)
      return
   end function
end program
Back to top
View user's profile Send private message
Display posts from previous:   
Post new topic   Reply to topic    forums.silverfrost.com Forum Index -> Support All times are GMT + 1 Hour
Page 1 of 1

 
Jump to:  
You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot vote in polls in this forum


Powered by phpBB © 2001, 2005 phpBB Group