Skip to content

Assumed shape vs contiguous arrays #46

@PierUgit

Description

@PierUgit

We already discussed that on Discourse.

BLAS/LAPACK routines mostly use assumed size arrays, thus contiguous arrays, but can handle non contiguous data by using the inc*/ld* arguments.

In the modern wrapper (e.g. solve) you are developing, you are defining assumed-shape arguments, which is a good thing.

The problem is that compilers may/will generate copy-in/copy-out in the underlying calls to the LAPACK routines.

subroutine brand_new(A)
   real, intent(inout) :: A(:,:)
   interface
      subroutine old_timer(n,A,lda)
         integer n, lda
         real A(lda,*)
      end subroutine old_timer
   end interface
   n = size(A,1)
   call old_timer(n,A,n)   ! possible copy-in/copy-out at this point
end subroutine brand_new

How can we avoid that? We talked about inquiring the strides, using the C interoperability and the C array descriptors. So it could look like:

subroutine brand_new(A)
   real, intent(inout) :: A(:,:)
   interface
      subroutine old_timer(n,A,lda)
         integer n, lda
         real A(lda,*)
      end subroutine old_timer
   end interface
   n = size(A,1)
   inc = stride(A,1)   ! C routine behind
   lda = stride(A,2)   ! C routine behind
   if (inc == 1) then
      call old_timer(n,A(1,1),lda)   ! not nice, but supposed to work and avoid copy-in/copy-out
   else
      call old_timer(n,A,n)   ! copy-in/copy-out cannot be avoided in this case
   end if
end subroutine brand_new

What do you think ?

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions