template_lapack_getrs.h

Go to the documentation of this file.
00001 /* Ergo, version 3.2, a program for linear scaling electronic structure
00002  * calculations.
00003  * Copyright (C) 2012 Elias Rudberg, Emanuel H. Rubensson, and Pawel Salek.
00004  * 
00005  * This program is free software: you can redistribute it and/or modify
00006  * it under the terms of the GNU General Public License as published by
00007  * the Free Software Foundation, either version 3 of the License, or
00008  * (at your option) any later version.
00009  * 
00010  * This program is distributed in the hope that it will be useful,
00011  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00012  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00013  * GNU General Public License for more details.
00014  * 
00015  * You should have received a copy of the GNU General Public License
00016  * along with this program.  If not, see <http://www.gnu.org/licenses/>.
00017  * 
00018  * Primary academic reference:
00019  * Kohn−Sham Density Functional Theory Electronic Structure Calculations 
00020  * with Linearly Scaling Computational Time and Memory Usage,
00021  * Elias Rudberg, Emanuel H. Rubensson, and Pawel Salek,
00022  * J. Chem. Theory Comput. 7, 340 (2011),
00023  * <http://dx.doi.org/10.1021/ct100611z>
00024  * 
00025  * For further information about Ergo, see <http://www.ergoscf.org>.
00026  */
00027  
00028  /* This file belongs to the template_lapack part of the Ergo source 
00029   * code. The source files in the template_lapack directory are modified
00030   * versions of files originally distributed as CLAPACK, see the
00031   * Copyright/license notice in the file template_lapack/COPYING.
00032   */
00033  
00034 
00035 #ifndef TEMPLATE_LAPACK_GETRS_HEADER
00036 #define TEMPLATE_LAPACK_GETRS_HEADER
00037 
00038 
00039 template<class Treal>
00040 int template_lapack_getrs(const char *trans, const integer *n, const integer *nrhs, 
00041         const Treal *a, const integer *lda, const integer *ipiv, Treal *b, const integer *
00042         ldb, integer *info)
00043 {
00044 /*  -- LAPACK routine (version 3.0) --   
00045        Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,   
00046        Courant Institute, Argonne National Lab, and Rice University   
00047        March 31, 1993   
00048 
00049 
00050     Purpose   
00051     =======   
00052 
00053     DGETRS solves a system of linear equations   
00054        A * X = B  or  A' * X = B   
00055     with a general N-by-N matrix A using the LU factorization computed   
00056     by DGETRF.   
00057 
00058     Arguments   
00059     =========   
00060 
00061     TRANS   (input) CHARACTER*1   
00062             Specifies the form of the system of equations:   
00063             = 'N':  A * X = B  (No transpose)   
00064             = 'T':  A'* X = B  (Transpose)   
00065             = 'C':  A'* X = B  (Conjugate transpose = Transpose)   
00066 
00067     N       (input) INTEGER   
00068             The order of the matrix A.  N >= 0.   
00069 
00070     NRHS    (input) INTEGER   
00071             The number of right hand sides, i.e., the number of columns   
00072             of the matrix B.  NRHS >= 0.   
00073 
00074     A       (input) DOUBLE PRECISION array, dimension (LDA,N)   
00075             The factors L and U from the factorization A = P*L*U   
00076             as computed by DGETRF.   
00077 
00078     LDA     (input) INTEGER   
00079             The leading dimension of the array A.  LDA >= max(1,N).   
00080 
00081     IPIV    (input) INTEGER array, dimension (N)   
00082             The pivot indices from DGETRF; for 1<=i<=N, row i of the   
00083             matrix was interchanged with row IPIV(i).   
00084 
00085     B       (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)   
00086             On entry, the right hand side matrix B.   
00087             On exit, the solution matrix X.   
00088 
00089     LDB     (input) INTEGER   
00090             The leading dimension of the array B.  LDB >= max(1,N).   
00091 
00092     INFO    (output) INTEGER   
00093             = 0:  successful exit   
00094             < 0:  if INFO = -i, the i-th argument had an illegal value   
00095 
00096     =====================================================================   
00097 
00098 
00099        Test the input parameters.   
00100 
00101        Parameter adjustments */
00102     /* Table of constant values */
00103      integer c__1 = 1;
00104      Treal c_b12 = 1.;
00105      integer c_n1 = -1;
00106     
00107     /* System generated locals */
00108     integer a_dim1, a_offset, b_dim1, b_offset, i__1;
00109     /* Local variables */
00110      logical notran;
00111 
00112 
00113     a_dim1 = *lda;
00114     a_offset = 1 + a_dim1 * 1;
00115     a -= a_offset;
00116     --ipiv;
00117     b_dim1 = *ldb;
00118     b_offset = 1 + b_dim1 * 1;
00119     b -= b_offset;
00120 
00121     /* Function Body */
00122     *info = 0;
00123     notran = template_blas_lsame(trans, "N");
00124     if (! notran && ! template_blas_lsame(trans, "T") && ! template_blas_lsame(
00125             trans, "C")) {
00126         *info = -1;
00127     } else if (*n < 0) {
00128         *info = -2;
00129     } else if (*nrhs < 0) {
00130         *info = -3;
00131     } else if (*lda < maxMACRO(1,*n)) {
00132         *info = -5;
00133     } else if (*ldb < maxMACRO(1,*n)) {
00134         *info = -8;
00135     }
00136     if (*info != 0) {
00137         i__1 = -(*info);
00138         template_blas_erbla("GETRS ", &i__1);
00139         return 0;
00140     }
00141 
00142 /*     Quick return if possible */
00143 
00144     if (*n == 0 || *nrhs == 0) {
00145         return 0;
00146     }
00147 
00148     if (notran) {
00149 
00150 /*        Solve A * X = B.   
00151 
00152           Apply row interchanges to the right hand sides. */
00153 
00154         template_lapack_laswp(nrhs, &b[b_offset], ldb, &c__1, n, &ipiv[1], &c__1);
00155 
00156 /*        Solve L*X = B, overwriting B with X. */
00157 
00158         template_blas_trsm("Left", "Lower", "No transpose", "Unit", n, nrhs, &c_b12, &a[
00159                 a_offset], lda, &b[b_offset], ldb);
00160 
00161 /*        Solve U*X = B, overwriting B with X. */
00162 
00163         template_blas_trsm("Left", "Upper", "No transpose", "Non-unit", n, nrhs, &c_b12, &
00164                 a[a_offset], lda, &b[b_offset], ldb);
00165     } else {
00166 
00167 /*        Solve A' * X = B.   
00168 
00169           Solve U'*X = B, overwriting B with X. */
00170 
00171         template_blas_trsm("Left", "Upper", "Transpose", "Non-unit", n, nrhs, &c_b12, &a[
00172                 a_offset], lda, &b[b_offset], ldb);
00173 
00174 /*        Solve L'*X = B, overwriting B with X. */
00175 
00176         template_blas_trsm("Left", "Lower", "Transpose", "Unit", n, nrhs, &c_b12, &a[
00177                 a_offset], lda, &b[b_offset], ldb);
00178 
00179 /*        Apply row interchanges to the solution vectors. */
00180 
00181         template_lapack_laswp(nrhs, &b[b_offset], ldb, &c__1, n, &ipiv[1], &c_n1);
00182     }
00183 
00184     return 0;
00185 
00186 /*     End of DGETRS */
00187 
00188 } /* dgetrs_ */
00189 
00190 #endif

Generated on Mon Sep 17 14:30:40 2012 for ergo by  doxygen 1.4.7