1 /*
2   Compute the Frobenius norm of a matrix.
3   */
4 #include <stdlib.h>
5 #include <stdio.h>
6 #include <math.h>
7 #include "declarations.h"
8 
Fnorm(A)9 double Fnorm(A)
10      struct blockmatrix A;
11 {
12   int blk;
13   double nrm;
14   double temp;
15 
16   nrm=0;
17   for (blk=1; blk<=A.nblocks; blk++)
18     {
19       switch (A.blocks[blk].blockcategory)
20 	{
21 	case DIAG:
22 	  temp=norm2(A.blocks[blk].blocksize,A.blocks[blk].data.vec+1);
23 	  nrm += temp*temp;
24 	  break;
25 	case MATRIX:
26 	  temp=norm2(A.blocks[blk].blocksize*A.blocks[blk].blocksize,
27 		     A.blocks[blk].data.mat);
28 	  nrm += temp*temp;
29 	  break;
30 	case PACKEDMATRIX:
31 	default:
32 	  printf("Fnorm illegal block type \n");
33 	  exit(206);
34 	};
35     };
36 
37   nrm=sqrt(nrm);
38   return(nrm);
39 }
40 
41 /*
42  * The Knorm is the sum of the Fnorms of the blocks.
43  */
44 
Knorm(A)45 double Knorm(A)
46      struct blockmatrix A;
47 {
48   int blk;
49   double nrm;
50   double temp;
51 
52   nrm=0;
53   for (blk=1; blk<=A.nblocks; blk++)
54     {
55       switch (A.blocks[blk].blockcategory)
56 	{
57 	case DIAG:
58 	  temp=norm2(A.blocks[blk].blocksize,A.blocks[blk].data.vec+1);
59 	  nrm += temp;
60 	  break;
61 	case MATRIX:
62 	  temp=norm2(A.blocks[blk].blocksize*A.blocks[blk].blocksize,
63 		     A.blocks[blk].data.mat);
64 	  nrm += temp;
65 	  break;
66 	case PACKEDMATRIX:
67 	default:
68 	  printf("Fnorm illegal block type \n");
69 	  exit(206);
70 	};
71     };
72 
73   return(nrm);
74 }
75 
76 
mat1norm(A)77 double mat1norm(A)
78      struct blockmatrix A;
79 {
80   int blk;
81   double nrm;
82   double temp;
83 
84   nrm=0;
85   for (blk=1; blk<=A.nblocks; blk++)
86     {
87       switch (A.blocks[blk].blockcategory)
88 	{
89 	case DIAG:
90 	  temp=norm1(A.blocks[blk].blocksize,A.blocks[blk].data.vec+1);
91 	  nrm += temp;
92 	  break;
93 	case MATRIX:
94 	  temp=norm1(1*A.blocks[blk].blocksize*A.blocks[blk].blocksize,
95 		     A.blocks[blk].data.mat);
96 	  nrm += temp;
97 	  break;
98 	case PACKEDMATRIX:
99 	default:
100 	  printf("Fnorm illegal block type \n");
101 	  exit(206);
102 	};
103     };
104 
105   return(nrm);
106 }
107 
108 
matinfnorm(A)109 double matinfnorm(A)
110      struct blockmatrix A;
111 {
112   int blk;
113   int i;
114   double nrm;
115 
116   nrm=0;
117   for (blk=1; blk<=A.nblocks; blk++)
118     {
119       switch (A.blocks[blk].blockcategory)
120 	{
121 	case DIAG:
122 	  for (i=1; i<=A.blocks[blk].blocksize; i++)
123 	    {
124 	      if (fabs(A.blocks[blk].data.vec[i]) > nrm)
125 		nrm=fabs(A.blocks[blk].data.vec[i]);
126 	    };
127 	  break;
128 	case MATRIX:
129 	  for (i=0; i<A.blocks[blk].blocksize*A.blocks[blk].blocksize; i++)
130 	    {
131 	      if (fabs(A.blocks[blk].data.vec[i]) > nrm)
132 		nrm=fabs(A.blocks[blk].data.vec[i]);
133 	    };
134 	  break;
135 	case PACKEDMATRIX:
136 	default:
137 	  printf("Fnorm illegal block type \n");
138 	  exit(206);
139 	};
140     };
141 
142   return(nrm);
143 }
144 
145