1
votes

I have a question about Fortran and correct allocation of allocatable user derived types.

Here is my code:

module polynom_mod
 implicit none

 type monomial
  integer,dimension(2) :: exponent
 end type

type polynom
  real, allocatable, dimension(:) :: coeff
  type(monomial),allocatable, dimension(:)   :: monom
  logical :: allocated
 !recursive type
  type(polynom),pointer :: p_dx,p_dy
 contains
  procedure :: init
  procedure :: init_dx
end type

here I want to derive a type polynom where I can do things like:

p%coeff(1)=1.0 
p%monom(1)%exponent(1)=2

and something like:

p%p_dx%coeff(1)=1.0 
p%p_dx%monom(1)%exponent(1)=2

so I wrote some init type-bound procedures where I can initialize my and allocate my types:

contains

function init(this,num) result(stat)
  implicit none
  integer, intent(in)      :: num
  class(polynom),intent(inout) :: this
  logical :: stat

  allocate(this%coeff(num))
  allocate(this%monom(num))

  this%allocated = .TRUE.
  stat = .TRUE.
end function

function init_dx(this,num) result(stat)
  implicit none

  integer, intent(in)      :: num
  class(polynom),intent(inout) :: this

  logical :: stat

  allocate(this%p_dx%coeff(num))
  allocate(this%p_dx%monom(num))

  this%p_dx%allocated = .TRUE.
  stat = .TRUE.
 end function   
end module

program testpolytype
 use polynom_mod

 type(polynom) :: p

 if(p%init(2)) then
  print *,"Polynom allocated!"
 end if

 if(p%p_dx%init_dx(2)) then
  print *,"Polynom_dx allocated!"
 end if

end program

This will compile with gfortran 4.6.3 but when I ran it I got a segmentation fault!

Is there a way to allocate recursive allocatable types?

2

2 Answers

2
votes

The superficial problem with your code is that, when the expression p%p_dx%init_dx(2) is computed the pointer component p%p_dx is undefined, and the segmentation fault is raised. Note that the pointer is undefined and not just not associated.

Right now I'm struggling to come up with a quick fix. The long fix would be to address what I think is a serious flaw in your approach; note that this is my opinion rather than a matter of black or white so read on only if you care for my input.

The functions init and init_dx are not free of side-effects, indeed they could be said to be almost all side-effect -- they return a logical value and, as a side-effect, initialise a polynom variable. The program seems to have no way to initialise a polynom without evaluating init and no way to evaluate init without wrapping it into a statement such as

if (p%init(2)) then
end if

You could, I suppose, rewrite these initialisation functions as subroutines, perhaps with a signature such as

call initialise_polynom(p,2)

This would, at least, remove the stain of impure functions from your code. But a better approach would be to write a function such as:

function new_poly(num)
  implicit none
  integer, intent(in) :: num
  type(polynom) :: new_poly
  allocate(new_poly%coeff(num))
  allocate(new_poly%monom(num))
  allocate(new_poly%p_dx)
end function new_poly

which

a) returns a new polynom; and

b) allocates component p_dx; and

c) is side-effect free.

You can then create a new polynom with an expression such as

p = new_poly(3)

and initialise the component with an expression such as

p%p_dx = new_poly(3)
0
votes

Answering my own question, I came up with an other solution witch also works without pointers but it's not so elegant like Marks one.

Define an other type :

type p_dx
 real, allocatable, dimension(:) :: coeff
 type(monomial),allocatable, dimension(:)   :: monom
 logical :: allocated
end type

and then use this with :

type polynom
 real, allocatable, dimension(:) :: coeff
 type(monomial),allocatable, dimension(:)   :: monom
 type(p_dx) :: dx
 logical :: allocated
contains
 procedure     :: init
end type

so you can do something like:

type(polynom) :: p

p%init(2)
p%dx%init_dx(3)