From dc2695a11f8a15bc83a8515a62ca6c234101c9f3 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?S=C3=A9bastien=20Villemot?= <sebastien@dynare.org>
Date: Thu, 3 Jun 2021 19:00:09 +0200
Subject: [PATCH] mjdgges and block_trust_region MEX: optimise by marking some
 function arguments as contiguous

This avoids unnecessary array copies before calling BLAS/LAPACK functions.
---
 mex/sources/block_trust_region/trust_region.f08 | 11 +++++++----
 mex/sources/mjdgges/mjdgges.F08                 |  6 ++++--
 2 files changed, 11 insertions(+), 6 deletions(-)

diff --git a/mex/sources/block_trust_region/trust_region.f08 b/mex/sources/block_trust_region/trust_region.f08
index 84e0f1c6cb..60e737d3cf 100644
--- a/mex/sources/block_trust_region/trust_region.f08
+++ b/mex/sources/block_trust_region/trust_region.f08
@@ -2,7 +2,7 @@
 !
 ! Implementation heavily inspired from the hybrj function from MINPACK
 
-! Copyright © 2019 Dynare Team
+! Copyright © 2019-2021 Dynare Team
 !
 ! This file is part of Dynare.
 !
@@ -223,11 +223,14 @@ contains
   ! designates element-by-element multiplication),
   ! x is a convex combination of the Gauss-Newton and scaled gradient
   subroutine dogleg(r, b, d, delta, x, gn, recompute_gn)
-    real(real64), dimension(:), intent(in) :: b, d
-    real(real64), dimension(:,:), intent(in) :: r
+    ! The arrays used in BLAS/LAPACK calls are required to be contiguous, to
+    ! avoid temporary copies before calling BLAS/LAPACK.
+    real(real64), dimension(:), contiguous, intent(in) :: b
+    real(real64), dimension(:), intent(in) :: d
+    real(real64), dimension(:,:), contiguous, intent(in) :: r
     real(real64), intent(in) :: delta ! Radius of the trust region
     real(real64), dimension(:), intent(out) :: x ! Solution of the problem
-    real(real64), dimension(:), intent(inout) :: gn ! Gauss-Newton direction
+    real(real64), dimension(:), contiguous, intent(inout) :: gn ! Gauss-Newton direction
     logical, intent(in) :: recompute_gn ! Whether to re-compute Gauss-Newton direction
 
     integer(blint) :: n
diff --git a/mex/sources/mjdgges/mjdgges.F08 b/mex/sources/mjdgges/mjdgges.F08
index 75c3f658dc..74de12568f 100644
--- a/mex/sources/mjdgges/mjdgges.F08
+++ b/mex/sources/mjdgges/mjdgges.F08
@@ -18,7 +18,7 @@
 !   eigval       [complex] (n×1) vector of generalized eigenvalues
 !   info         [integer] scalar, error code of dgges (or 30 if eigenvalue close to 0÷0)
 
-! Copyright © 2006-2020 Dynare Team
+! Copyright © 2006-2021 Dynare Team
 !
 ! This file is part of Dynare.
 !
@@ -69,7 +69,9 @@ subroutine mexFunction(nlhs, plhs, nrhs, prhs) bind(c, name='mexFunction')
   integer(blint) :: n_bl, lwork, info_bl, sdim_bl
   real(real64), dimension(:), allocatable :: alpha_r, alpha_i, beta, work
   logical(bllog), dimension(:), allocatable :: bwork
-  real(real64), dimension(:), pointer :: s, t, z, info, sdim, vsl
+  ! The pointers used in the LAPACK call are marked as contiguous, to
+  ! avoid temporary copies beforehand.
+  real(real64), dimension(:), pointer, contiguous :: s, t, z, info, sdim, vsl
 #if MX_HAS_INTERLEAVED_COMPLEX
   complex(real64), dimension(:), pointer :: gev
 #else
-- 
GitLab