1/* 2 luserver.h 3 4 LU factorization Web service. 5 6-------------------------------------------------------------------------------- 7gSOAP XML Web services tools 8Copyright (C) 2001-2008, Robert van Engelen, Genivia, Inc. All Rights Reserved. 9This software is released under one of the following two licenses: 10GPL or Genivia's license for commercial use. 11-------------------------------------------------------------------------------- 12GPL license. 13 14This program is free software; you can redistribute it and/or modify it under 15the terms of the GNU General Public License as published by the Free Software 16Foundation; either version 2 of the License, or (at your option) any later 17version. 18 19This program is distributed in the hope that it will be useful, but WITHOUT ANY 20WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A 21PARTICULAR PURPOSE. See the GNU General Public License for more details. 22 23You should have received a copy of the GNU General Public License along with 24this program; if not, write to the Free Software Foundation, Inc., 59 Temple 25Place, Suite 330, Boston, MA 02111-1307 USA 26 27Author contact information: 28engelen@genivia.com / engelen@acm.org 29-------------------------------------------------------------------------------- 30A commercial use license is available from Genivia, Inc., contact@genivia.com 31-------------------------------------------------------------------------------- 32*/ 33 34#include "soapH.h" 35#include <math.h> 36 37//////////////////////////////////////////////////////////////////////////////// 38// 39// main 40// 41//////////////////////////////////////////////////////////////////////////////// 42 43int main(int argc, char **argv) 44{ struct soap *soap; 45 int m, s; 46 soap = soap_new(); 47 if (argc < 3) 48 soap_serve(soap); // run as CGI application over the Web 49 else // run as stand-alone server on machine given by argv[1] listening to port argv[2] 50 { m = soap_bind(soap, argv[1], atoi(argv[2]), 100); 51 if (m < 0) 52 { soap_print_fault(soap, stderr); 53 exit(-1); 54 } 55 fprintf(stderr, "Socket connection successful: master socket = %d\n", m); 56 for (int i = 1; ; i++) 57 { s = soap_accept(soap); 58 if (s < 0) 59 { soap_print_fault(soap, stderr); 60 exit(-1); 61 } 62 fprintf(stderr, "%d: accepted connection from IP = %d.%d.%d.%d socket = %d ... ", i, (int)(soap->ip>>24)&0xFF, (int)(soap->ip>>16)&0xFF, (int)(soap->ip>>8)&0xFF, (int)soap->ip&0xFF, s); 63 soap_serve(soap); // process request 64 fprintf(stderr, "request served\n"); 65 soap_destroy(soap); // delete class instances 66 soap_end(soap); // clean up everything and close socket 67 } 68 } 69 soap_done(soap); 70 return 0; 71} 72 73//////////////////////////////////////////////////////////////////////////////// 74// 75// LU decomposition: remote method interface 76// 77//////////////////////////////////////////////////////////////////////////////// 78 79int ludcmp(struct soap*, matrix&, ivector&, double&); 80 81int ns1__ludcmp(struct soap *soap, matrix *a, struct ns1__ludcmpResponse &result) 82{ result.a = a; 83 result.i = soap_new_ivector(soap, -1); 84 if (ludcmp(soap, *result.a, *result.i, result.d)) 85 return soap_sender_fault(soap, "Singular matrix in routine ludcmp", NULL); 86 return SOAP_OK; 87} 88 89//////////////////////////////////////////////////////////////////////////////// 90// 91// LU decomposition: algorithm 92// 93//////////////////////////////////////////////////////////////////////////////// 94 95int ludcmp(struct soap *soap, matrix &a, ivector &indx, double &d) 96{ int i, imax = 0, j, k, n; 97 double big, dum, sum, temp; 98 n = a.size(); 99 vector vv(soap); 100 vv.resize(n); 101 indx.resize(n); 102 d = 1.0; 103 for (i = 1; i <= n; i++) 104 { big = 0.0; 105 a[i].resize(n); 106 for (j = 1; j <= n; j++) 107 if ((temp = fabs(a[i][j])) > big) 108 big = temp; 109 if (big == 0.0) 110 return -1; 111 vv[i] = 1.0/big; 112 } 113 for (j = 1; j <= n; j++) 114 { for (i = 1; i < j; i++) 115 { sum = a[i][j]; 116 for (k = 1; k < i; k++) 117 sum -= a[i][k]*a[k][j]; 118 a[i][j] = sum; 119 } 120 big = 0.0; 121 for (i = j; i <= n; i++) 122 { sum = a[i][j]; 123 for (k = 1; k < j; k++) 124 sum -= a[i][k]*a[k][j]; 125 a[i][j] = sum; 126 if ((dum = vv[i]*fabs(sum)) >= big) 127 { big = dum; 128 imax = i; 129 } 130 } 131 if (j != imax) 132 { for (k = 1; k <= n; k++) 133 { dum = a[imax][k]; 134 a[imax][k] = a[j][k]; 135 a[j][k] = dum; 136 } 137 d = -d; 138 vv[imax] = vv[j]; 139 } 140 indx[j] = imax; 141 if (a[j][j] == 0.0) 142 a[j][j] = 1.0e-20; 143 if (j != n) 144 { dum = 1.0/a[j][j]; 145 for (i = j+1; i <= n; i++) 146 a[i][j] *= dum; 147 } 148 } 149 for (i = 1; i <= n; i++) 150 { for (j = 1; j <= n; j++) 151 if (fabs(a[i][j]) > 1.0e-15) 152 break; 153 for (k = n; k > j; k--) 154 if (fabs(a[i][k]) > 1.0e-15) 155 break; 156 a[i].resize(j, k); 157 } 158 return 0; 159} 160 161//////////////////////////////////////////////////////////////////////////////// 162// 163// Forward- and backsubstitution: remote method interface 164// 165//////////////////////////////////////////////////////////////////////////////// 166 167int lubksb(matrix&, ivector&, vector &b); 168 169int ns1__lubksb(struct soap *soap, matrix *a, ivector *i, vector *b, vector *x) 170{ lubksb(*a, *i, *b); 171 *x = *b; 172 return SOAP_OK; 173} 174 175//////////////////////////////////////////////////////////////////////////////// 176// 177// Forward- and backsubstitution: algorithm 178// 179//////////////////////////////////////////////////////////////////////////////// 180 181int lubksb(matrix &a, ivector &indx, vector &b) 182{ int i, j, k, ip, n, m, ii = 0; 183 double sum; 184 n = a.size(); 185 b.resize(n); 186 for (i = 1; i <= n; i++) 187 { ip = indx[i]; 188 sum = b[ip]; 189 b[ip] = b[i]; 190 if (ii) 191 { k = a[i].start(); 192 if (ii > k) 193 k = ii; 194 m = a[i].end(); 195 if (i-1 < m) 196 m = i-1; 197 for (j = k; j <= m; j++) 198 sum -= a[i][j]*b[j]; 199 } 200 else if (sum) 201 ii = i; 202 b[i] = sum; 203 } 204 for (i = n; i >= 1; i--) 205 { sum = b[i]; 206 k = a[i].start(); 207 if (i+1 > k) 208 k = i+1; 209 m = a[i].end(); 210 if (n < m) 211 m = n; 212 for (j = k; j <= m; j++) 213 sum -= a[i][j]*b[j]; 214 b[i] = sum/a[i][i]; 215 } 216 return SOAP_OK; 217} 218 219//////////////////////////////////////////////////////////////////////////////// 220// 221// LU solve: remote method interface 222// 223//////////////////////////////////////////////////////////////////////////////// 224 225int ns1__lusol(struct soap *soap, matrix *a, vector *b, vector *x) 226{ ivector indx(soap); 227 double d; 228 if (ludcmp(soap, *a, indx, d)) 229 return soap_sender_fault(soap, "Singular matrix in routine ludcmp", NULL); 230 lubksb(*a, indx, *b); 231 *x = *b; 232 return SOAP_OK; 233} 234 235//////////////////////////////////////////////////////////////////////////////// 236// 237// LU solve multiple: remote method interface 238// 239//////////////////////////////////////////////////////////////////////////////// 240 241int ns1__lusols(struct soap *soap, matrix *a, matrix *b, matrix *x) 242{ ivector indx(soap); 243 double d; 244 if (ludcmp(soap, *a, indx, d)) 245 return soap_sender_fault(soap, "Singular matrix in routine ludcmp", NULL); 246 for (int i = 1; i <= b->size(); i++) 247 lubksb(*a, indx, (*b)[i]); 248 *x = *b; 249 return SOAP_OK; 250} 251 252//////////////////////////////////////////////////////////////////////////////// 253// 254// LU inverse: remote method interface 255// 256//////////////////////////////////////////////////////////////////////////////// 257 258int ns1__luinv(struct soap *soap, matrix *a, matrix *b) 259{ vector col(soap); 260 ivector indx(soap); 261 double d; 262 int i, j, k, n; 263 if (ludcmp(soap, *a, indx, d)) 264 return soap_sender_fault(soap, "Singular matrix in routine ludcmp", NULL); 265 n = a->size(); 266 col.resize(n); 267 b->resize(n, n); 268 for (j = 1; j <= n; j++) 269 { for (i = 1; i <= n; i++) 270 col[i] = 0.0; 271 col[j] = 1.0; 272 lubksb(*a, indx, col); 273 for (i = 1; i <= n; i++) 274 (*b)[i][j] = col[i]; 275 } 276 for (i = 1; i <= n; i++) 277 { for (j = 1; j <= n; j++) 278 if (fabs((*b)[i][j]) > 1.0e-15) 279 break; 280 for (k = n; k > j; k--) 281 if (fabs((*b)[i][k]) > 1.0e-15) 282 break; 283 (*b)[i].resize(j, k); 284 } 285 return SOAP_OK; 286} 287 288//////////////////////////////////////////////////////////////////////////////// 289// 290// LU determinant: remote method interface 291// 292//////////////////////////////////////////////////////////////////////////////// 293 294int ns1__ludet(struct soap *soap, matrix *a, double &d) 295{ ivector indx(soap); 296 if (ludcmp(soap, *a, indx, d)) 297 return soap_sender_fault(soap, "Singular matrix in routine ludcmp", NULL); 298 for (int i = 1; i <= a->__size; i++) 299 d *= (*a)[i][i]; 300 return SOAP_OK; 301} 302 303struct Namespace namespaces[] = 304{// "ns-prefix", "ns-name", "ns-pattern" 305 {"SOAP-ENV", "http://schemas.xmlsoap.org/soap/envelope/"}, 306 {"SOAP-ENC", "http://schemas.xmlsoap.org/soap/encoding/"}, 307 {"xsi", "http://www.w3.org/2001/XMLSchema-instance", "http://wwww.w3.org/*/XMLSchema-instance"}, 308 {"xsd", "http://www.w3.org/2001/XMLSchema", "http://www.w3.org/*/XMLSchema"}, 309 {"ns1", "urn:lu"}, 310 {NULL, NULL} 311}; 312