operators for matrix

Dear all,

I am writeing a program which will contain many matrix manipulations.
I want to define my own operators to shorten the syntax for the matrix
manipulations. Below is a small program as example.

This program run fine using g95.

Using the ifort 9.1 compiler give error:
Syntax error, found '.' when expecting one of: ( , <END-OF-STATEMENT> ;
<IDENTIFIER> <CHAR_CON_KIND_PARAM> <CHAR_NAM_KIND_PARAM> ...
write(*,*) .t.(a.mm..ms.(Im+(.ms.btot).mm.s)) After commenting out this
line it runs fine.

The gfortran compiler warns:
matrixmantst.f90: In function ‘vectorskew’:
matrixmantst.f90:55: warning: Function does not return a value
matrixmantst.f90: In function ‘matrixskew’:
matrixmantst.f90:47: warning: Function does not return a value
and gives a segmentation fault when it is run. If I changes:
call skewmatrix(MatrixSkew,vec1(:,1))
by:
real(dp) :: tmm(3,3)
call skewmatrix(tmm,vec1(:,1))
MatrixSkew = tmm
(same for VectorSkew)
it runs fine without runtime checks, with runtime check it gives still
a segmentation fault

The Sun f95 compiler gives error:
public operator(.mm.),operator(.t.),operator(.ms.), skewmatrix
^
"matrixmantst.f90", Line = 10, Column = 32: ERROR: A defined operator
must not be the same as a logical literal constant. After changing t by
tt it runs fine.

The PGI compiler gives:
PGF90-S-0034-Syntax error at or near ) (matrixmantst.f90: 102)
PGF90-S-0034-Syntax error at or near ) (matrixmantst.f90: 108)
By adding additional parentheses (.t.(a.mm.(.ms.(Im+(.ms.btot).mm.s))))
this error was is resolved and the program runs fine.


My questions are: Is this code according the F95 standards? Have you
any suggestions/comments how I can improve this code. Especially, how
can this be done so it will work for every F95 compiler?

Thanks in advance,

Johan



module types
! select kind
integer, parameter :: dp=kind(0.d0)
end module

module matrixman
! module for short matrix syntax
private

public operator(.mm.),operator(.t.),operator(.ms.), skewmatrix

interface operator (.mm.)
module procedure MatrixMultiplication
end interface

interface operator (.t.)
module procedure MatrixTranspose
end interface

interface operator (.ms.)
module procedure MatrixSkew
module procedure VectorSkew
end interface

contains
function MatrixMultiplication(mat1,mat2)
use types
implicit none
real(dp), intent(in) :: mat1(:,:), mat2(:,:)
real(dp) :: MatrixMultiplication(size(mat1,dim=1),size(mat2,dim=2))
MatrixMultiplication = matmul(mat1,mat2)
end function

function MatrixTranspose(mat1)
use types
implicit none
real(dp), intent(in) :: mat1(:,:)
real(dp) :: MatrixTranspose(size(mat1,dim=2),size(mat1,dim=1))
MatrixTranspose = transpose(mat1)
end function

function MatrixSkew(vec1)
use types
implicit none
real(dp), intent(in) :: vec1(3,1)
real(dp) :: MatrixSkew(3,3)
call skewmatrix(MatrixSkew,vec1(:,1))
end function

function VectorSkew(vec1)
use types
implicit none
real(dp), intent(in) :: vec1(3)
real(dp) :: VectorSkew(3,3)
call skewmatrix(VectorSkew,vec1)
end function

subroutine skewmatrix(SMat,s)
use types
implicit none

real(dp), intent(in) :: s(3)
real(dp), intent(out) :: SMat(3,3)

SMat(1,1) = 0.0_dp
SMat(1,2) = -s(3)
SMat(1,3) = s(2)

SMat(2,1) = s(3)
SMat(2,2) = 0.0_dp
SMat(2,3) = -s(1)

SMat(3,1) = -s(2)
SMat(3,2) = s(1)
SMat(3,3) = 0.0_dp
end subroutine
end module


program matrixtst
use types
use matrixman
implicit none
real(dp) :: res1(3,3), a(3,3), Im(3,3), btot(3,1), s(3,1)
real(dp) :: res2(3,3), bk(3,3), w1(3,1)

! filling the matrix with random numbers
call random_seed()
call random_number(a)
call random_number(Im)
call random_number(btot)
call random_number(s)

! original code
call skewmatrix(bk,btot)
bk = bk + Im
w1 = matmul(bk,s)
call skewmatrix(bk,w1)
res1 = transpose(matmul(a,bk))

! short code which should be equivalent
res2 = .t.(a.mm..ms.(Im+(.ms.btot).mm.s))

write(*,*) res1
write(*,*)
write(*,*) res2
write(*,*)
write(*,*) .t.(a.mm..ms.(Im+(.ms.btot).mm.s))

end program
.

Relevant Pages

  • Re: Ada exception block does NOT work?
    ... The parameter/formal declarations are like ... the compiler would guess for you how to call puts. ... > * implicit declaration of external functions. ... > K&R probably where good assember programmers but they where lousy compiler ...
    (comp.lang.ada)
  • Re: trim(string) problems
    ... >> Together these make me think it is a compiler problem ... >> but correctly finds size with both explicit and implicit ... that makes my test program with the module work and should ...
    (comp.lang.fortran)
  • Re: first time using direct access files: recl ??
    ... I suppose the compiler *MIGHT* be omitting non-written trailing ... he was particularly metioning source form issues. ... implicit typing (evil), the compiler will just make one up anyway. ...
    (comp.lang.fortran)
  • Re: C# 3.0 Proposals
    ... Implicit means that it is assumed what the type is, ... Whereas with var, you can't do it, you will get a compiler error. ... "26.1 Implicitly typed local variables ...
    (microsoft.public.dotnet.languages.csharp)