$search
00001 // Copyright (C) 2008-2011 NICTA (www.nicta.com.au) 00002 // Copyright (C) 2008-2011 Conrad Sanderson 00003 // 00004 // This file is part of the Armadillo C++ library. 00005 // It is provided without any warranty of fitness 00006 // for any purpose. You can redistribute this file 00007 // and/or modify it under the terms of the GNU 00008 // Lesser General Public License (LGPL) as published 00009 // by the Free Software Foundation, either version 3 00010 // of the License or (at your option) any later version. 00011 // (see http://www.opensource.org/licenses for more info) 00012 00013 00016 00017 00018 00019 template<typename T1> 00020 inline 00021 void 00022 op_diagmat::apply(Mat<typename T1::elem_type>& out, const Op<T1, op_diagmat>& X) 00023 { 00024 arma_extra_debug_sigprint(); 00025 00026 typedef typename T1::elem_type eT; 00027 00028 const unwrap<T1> tmp(X.m); 00029 const Mat<eT>& A = tmp.M; 00030 00031 if(A.is_vec() == true) 00032 { 00033 // generate a diagonal matrix out of a vector 00034 00035 const uword N = A.n_elem; 00036 const eT* A_mem = A.memptr(); 00037 00038 if(&out != &A) 00039 { 00040 // no aliasing 00041 out.zeros(N,N); 00042 00043 for(uword i=0; i<N; ++i) 00044 { 00045 out.at(i,i) = A_mem[i]; 00046 } 00047 } 00048 else 00049 { 00050 // aliasing 00051 00052 const podarray<eT> tmp(A_mem, N); 00053 00054 const eT* tmp_mem = tmp.memptr(); 00055 00056 out.zeros(N,N); 00057 00058 for(uword i=0; i<N; ++i) 00059 { 00060 out.at(i,i) = tmp_mem[i]; 00061 } 00062 } 00063 } 00064 else 00065 { 00066 // generate a diagonal matrix out of a matrix 00067 00068 arma_debug_check( (A.is_square() == false), "diagmat(): given matrix is not square" ); 00069 00070 const uword N = A.n_rows; 00071 00072 if(&out != &A) 00073 { 00074 // no aliasing 00075 00076 out.zeros(N,N); 00077 00078 for(uword i=0; i<N; ++i) 00079 { 00080 out.at(i,i) = A.at(i,i); 00081 } 00082 } 00083 else 00084 { 00085 // aliasing 00086 00087 for(uword i=0; i<N; ++i) 00088 { 00089 eT* colptr = out.colptr(i); 00090 00091 // clear above the diagonal 00092 arrayops::inplace_set(colptr, eT(0), i); 00093 00094 // clear below the diagonal 00095 arrayops::inplace_set(colptr+(i+1), eT(0), N-1-i); 00096 } 00097 } 00098 } 00099 } 00100 00101 00102