initializing diagonal matrix with some vector

[David Necas (Yeti) <yeti at physics dot muni dot cz>]

Hello,

I often need to initialize a diagonal matrix from
a vector.  I think a function for such an
initialization could be generally useful, so
here's a patch, if you think so

Yeti


--- gsl-1.1.1.orig/matrix/gsl_matrix_char.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_char.h	Sun Apr 28 02:35:55 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_char_set_zero (gsl_matrix_char * m);
 void gsl_matrix_char_set_identity (gsl_matrix_char * m);
+void gsl_matrix_char_set_diagonal (gsl_matrix_char * m, gsl_vector_char * v);
 void gsl_matrix_char_set_all (gsl_matrix_char * m, char x);
 
 int gsl_matrix_char_fread (FILE * stream, gsl_matrix_char * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_complex_double.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_complex_double.h	Sun Apr 28 02:36:12 2002
@@ -191,6 +191,7 @@
 
 void gsl_matrix_complex_set_zero (gsl_matrix_complex * m);
 void gsl_matrix_complex_set_identity (gsl_matrix_complex * m);
+void gsl_matrix_complex_set_diagonal (gsl_matrix_complex * m, gsl_vector_complex * v);
 void gsl_matrix_complex_set_all (gsl_matrix_complex * m, gsl_complex x);
 
 int gsl_matrix_complex_fread (FILE * stream, gsl_matrix_complex * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_complex_float.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_complex_float.h	Sun Apr 28 02:39:52 2002
@@ -191,6 +191,7 @@
 
 void gsl_matrix_complex_float_set_zero (gsl_matrix_complex_float * m);
 void gsl_matrix_complex_float_set_identity (gsl_matrix_complex_float * m);
+void gsl_matrix_complex_float_set_diagonal (gsl_matrix_complex_float * m, gsl_vector_complex_float * v);
 void gsl_matrix_complex_float_set_all (gsl_matrix_complex_float * m, gsl_complex_float x);
 
 int gsl_matrix_complex_float_fread (FILE * stream, gsl_matrix_complex_float * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_complex_long_double.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_complex_long_double.h	Sun Apr 28 02:39:27 2002
@@ -191,6 +191,7 @@
 
 void gsl_matrix_complex_long_double_set_zero (gsl_matrix_complex_long_double * m);
 void gsl_matrix_complex_long_double_set_identity (gsl_matrix_complex_long_double * m);
+void gsl_matrix_complex_long_double_set_diagonal (gsl_matrix_complex_long_double * m, gsl_vector_complex_long_double * v);
 void gsl_matrix_complex_long_double_set_all (gsl_matrix_complex_long_double * m, gsl_complex_long_double x);
 
 int gsl_matrix_complex_long_double_fread (FILE * stream, gsl_matrix_complex_long_double * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_double.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_double.h	Sun Apr 28 02:37:05 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_set_zero (gsl_matrix * m);
 void gsl_matrix_set_identity (gsl_matrix * m);
+void gsl_matrix_set_diagonal (gsl_matrix * m, gsl_vector * v);
 void gsl_matrix_set_all (gsl_matrix * m, double x);
 
 int gsl_matrix_fread (FILE * stream, gsl_matrix * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_float.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_float.h	Sun Apr 28 02:37:13 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_float_set_zero (gsl_matrix_float * m);
 void gsl_matrix_float_set_identity (gsl_matrix_float * m);
+void gsl_matrix_float_set_diagonal (gsl_matrix_float * m, gsl_vector_float * v);
 void gsl_matrix_float_set_all (gsl_matrix_float * m, float x);
 
 int gsl_matrix_float_fread (FILE * stream, gsl_matrix_float * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_int.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_int.h	Sun Apr 28 02:37:18 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_int_set_zero (gsl_matrix_int * m);
 void gsl_matrix_int_set_identity (gsl_matrix_int * m);
+void gsl_matrix_int_set_diagonal (gsl_matrix_int * m, gsl_vector_int * v);
 void gsl_matrix_int_set_all (gsl_matrix_int * m, int x);
 
 int gsl_matrix_int_fread (FILE * stream, gsl_matrix_int * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_long_double.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_long_double.h	Sun Apr 28 02:40:14 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_long_double_set_zero (gsl_matrix_long_double * m);
 void gsl_matrix_long_double_set_identity (gsl_matrix_long_double * m);
+void gsl_matrix_long_double_set_diagonal (gsl_matrix_long_double * m, gsl_vector_long_double * v);
 void gsl_matrix_long_double_set_all (gsl_matrix_long_double * m, long double x);
 
 int gsl_matrix_long_double_fread (FILE * stream, gsl_matrix_long_double * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_long.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_long.h	Sun Apr 28 02:37:27 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_long_set_zero (gsl_matrix_long * m);
 void gsl_matrix_long_set_identity (gsl_matrix_long * m);
+void gsl_matrix_long_set_diagonal (gsl_matrix_long * m, gsl_vector_long * v);
 void gsl_matrix_long_set_all (gsl_matrix_long * m, long x);
 
 int gsl_matrix_long_fread (FILE * stream, gsl_matrix_long * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_short.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_short.h	Sun Apr 28 02:37:31 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_short_set_zero (gsl_matrix_short * m);
 void gsl_matrix_short_set_identity (gsl_matrix_short * m);
+void gsl_matrix_short_set_diagonal (gsl_matrix_short * m, gsl_vector_short * v);
 void gsl_matrix_short_set_all (gsl_matrix_short * m, short x);
 
 int gsl_matrix_short_fread (FILE * stream, gsl_matrix_short * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_uchar.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_uchar.h	Sun Apr 28 02:37:35 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_uchar_set_zero (gsl_matrix_uchar * m);
 void gsl_matrix_uchar_set_identity (gsl_matrix_uchar * m);
+void gsl_matrix_uchar_set_diagonal (gsl_matrix_uchar * m, gsl_vector_uchar * v);
 void gsl_matrix_uchar_set_all (gsl_matrix_uchar * m, unsigned char x);
 
 int gsl_matrix_uchar_fread (FILE * stream, gsl_matrix_uchar * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_uint.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_uint.h	Sun Apr 28 02:37:40 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_uint_set_zero (gsl_matrix_uint * m);
 void gsl_matrix_uint_set_identity (gsl_matrix_uint * m);
+void gsl_matrix_uint_set_diagonal (gsl_matrix_uint * m, gsl_vector_uint * v);
 void gsl_matrix_uint_set_all (gsl_matrix_uint * m, unsigned int x);
 
 int gsl_matrix_uint_fread (FILE * stream, gsl_matrix_uint * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_ulong.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_ulong.h	Sun Apr 28 02:37:46 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_ulong_set_zero (gsl_matrix_ulong * m);
 void gsl_matrix_ulong_set_identity (gsl_matrix_ulong * m);
+void gsl_matrix_ulong_set_diagonal (gsl_matrix_ulong * m, gsl_vector_ulong * v);
 void gsl_matrix_ulong_set_all (gsl_matrix_ulong * m, unsigned long x);
 
 int gsl_matrix_ulong_fread (FILE * stream, gsl_matrix_ulong * m) ;
--- gsl-1.1.1.orig/matrix/gsl_matrix_ushort.h	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/gsl_matrix_ushort.h	Sun Apr 28 02:37:55 2002
@@ -194,6 +194,7 @@
 
 void gsl_matrix_ushort_set_zero (gsl_matrix_ushort * m);
 void gsl_matrix_ushort_set_identity (gsl_matrix_ushort * m);
+void gsl_matrix_ushort_set_diagonal (gsl_matrix_ushort * m, gsl_vector_ushort * v);
 void gsl_matrix_ushort_set_all (gsl_matrix_ushort * m, unsigned short x);
 
 int gsl_matrix_ushort_fread (FILE * stream, gsl_matrix_ushort * m) ;
--- gsl-1.1.1.orig/matrix/init_source.c	Wed May  1 15:42:07 2002
+++ gsl-1.1.1/matrix/init_source.c	Sun Apr 28 02:31:02 2002
@@ -211,6 +211,25 @@
 }
 
 void
+FUNCTION (gsl_matrix, set_diagonal) (TYPE (gsl_matrix) * m, TYPE (gsl_vector) * v)
+{
+  size_t i, j;
+  ATOMIC * const data = m->data;
+  const size_t p = m->size1 ;
+  const size_t q = m->size2 ;
+  const size_t tda = m->tda ;
+
+  for (i = 0; i < p; i++)
+    {
+      for (j = 0; j < q; j++)
+        {
+          *(BASE *) (data + MULTIPLICITY * (i * tda + j))
+            = ((i == j) ? *(BASE *) (v->data + MULTIPLICITY * i * v->stride) : zero);
+        }
+    }
+}
+
+void
 FUNCTION (gsl_matrix, set_zero) (TYPE (gsl_matrix) * m)
 {
   size_t i, j;