📄 dsrc2c.c
字号:
/** dsrc2c: THIS FILE HAS BEEN MODIFIED AFTER f2c
* Modifications were
* - removed zbrent - copyright issue
* - removed eqrt1s - copyright issue
* - In eigvns_, replaced call to eqrt1s with call to v3p_netlib_tqlrat_
* - removed eigvss_ - required zbrent
* - use eigvns_ instead
**/
#define V3P_NETLIB_SRC
#include "v3p_netlib.h"
/* Modified by Peter Vanroose, Oct 2003: manual optimisation and clean-up */
extern double log(double), sqrt(double); /* #include <math.h> */
extern long time(long *timer); /* #include <time.h> */
extern doublereal cheby_(doublereal *, doublereal *, doublereal *, integer *, doublereal *, doublereal *);
extern doublereal determ_(integer *, doublereal *, doublereal *);
extern doublereal eigvns_(integer *, doublereal *, doublereal *, doublereal *, integer *);
extern doublereal itpackddot_(integer *, doublereal *, integer *, doublereal *, integer *);
extern doublereal pbeta_(integer *, integer *, integer *, doublereal *, doublereal *, doublereal *, doublereal *);
extern doublereal pvtbv_(integer *, integer *, integer *, doublereal *, doublereal *);
extern doublereal tau_(integer *);
extern doublereal dsrc_timer_(real *);
extern integer bisrch_(integer *, integer *, integer *);
extern integer ipstr_(doublereal *);
extern logical chgsme_(doublereal *, integer *);
extern logical omgchg_(integer *);
extern logical omgstr_(integer *);
extern logical tstchg_(integer *);
/***** BEGIN VXL ADDITIONS ****/
/* Turn off warnings in f2c generated code */
#if defined(_MSC_VER)
# if defined(__ICL)
# pragma warning(disable: 239 264 1011 )
# else
# pragma warning(disable: 4101 4244 4554 4756 4723)
# endif
#endif
/***** END VXL ADDITIONS ****/
int echout_(integer *iparm, doublereal *rparm, integer *imthd);
int echall_(integer *nn, integer *ia, integer *ja, doublereal *a, doublereal *rhs,
integer *iparm, doublereal *rparm, integer *icall);
int ivfill_(integer *n, integer *iv, integer *ival);
int vfill_(integer *n, doublereal *v, doublereal *val);
int sbelm_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs, integer *iw,
doublereal *rw, doublereal *tol, integer *isym, integer *level, integer *nout, integer *ier);
int prbndx_(integer *n, integer *nblack, integer *ia, integer *ja, integer *p,
integer *ip, integer *level, integer *nout, integer *ier);
int permat_(integer *n, integer *ia, integer *ja, doublereal *a, integer *p,
integer *newia, integer *isym, integer * level, integer *nout, integer* ierr);
int pervec_(integer *n, doublereal *v, integer *p);
int scal_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs,
doublereal *u, doublereal *d, integer *level, integer *nout, integer *ier);
int itpackdcopy_(integer *n, doublereal *dx, integer *incx, doublereal *dy, integer *incy);
int pjac_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *u, doublereal *rhs);
int vevmw_(integer *n, doublereal *v, doublereal *w);
int perror_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs,
doublereal *u, doublereal *w, doublereal *digtt1, doublereal *digtt2, integer *idgtts);
int itjsi_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs,
doublereal *u, doublereal *u1, doublereal *d, integer *icnt);
int itsor_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs, doublereal *u, doublereal *wk);
int omeg_(doublereal *dnrm, integer *iflag);
int pfsor_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *u, doublereal *rhs);
int itsrcg_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs,
doublereal *u, doublereal *u1, doublereal *c, doublereal *c1,
doublereal *d, doublereal *dl, doublereal *wk, doublereal *tri);
int itsrsi_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs, doublereal *u,
doublereal *u1, doublereal *c, doublereal *d, doublereal *ctwd, doublereal *wk);
int chgcon_(doublereal *tri, doublereal *gamold, doublereal *rhoold, integer *ibmth);
int pstop_(integer *n, doublereal *u, doublereal *dnrm, doublereal *ccon, integer *iflag, logical *q1);
int parcon_(doublereal *dtnrm, doublereal *c1, doublereal *c2, doublereal *c3, doublereal *c4,
doublereal *gamold, doublereal * rhotmp, integer *ibmth);
int sum3_(integer *n, doublereal *c1, doublereal *x1, doublereal *c2, doublereal *x2, doublereal *c3, doublereal *x3);
int iterm_(integer *nn, doublereal *a, doublereal *u, doublereal *wk, integer *imthdd);
int chgsi_(doublereal *dtnrm, integer *ibmth);
int itpackdaxpy_(integer *n, doublereal *da, doublereal *dx, integer *incx, doublereal *dy, integer *incy);
int pfsor1_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *u, doublereal *rhs);
int vevpw_(integer *n, doublereal *v, doublereal *w);
int pbsor_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *u, doublereal *rhs);
int wevmw_(integer *n, doublereal *v, doublereal *w);
int pssor1_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *u,
doublereal *rhs, doublereal *fr, doublereal *br);
int pmult_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *u, doublereal *w);
int vout_(integer *n, doublereal *v, integer *iswt, integer *nout);
int unscal_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs, doublereal *u, doublereal *d);
int prsred_(integer *nb, integer *nr, integer *ia, integer *ja, doublereal *a, doublereal *ub, doublereal *vr);
int itrscg_(integer *n, integer *nb, integer *ia, integer *ja, doublereal *a, doublereal *ub,
doublereal *ub1, doublereal *db, doublereal *db1, doublereal *wb, doublereal *tri);
int itrssi_(integer *n, integer *nb, integer *ia, integer *ja, doublereal *a,
doublereal *rhs, doublereal *ub, doublereal *ub1, doublereal *db);
int parsi_(doublereal *c1, doublereal *c2, doublereal *c3, integer *ibmth);
int itjcg_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *u, doublereal *u1,
doublereal *d, doublereal *d1, doublereal *dtwd, doublereal *tri);
int prsblk_(integer *nb, integer *nr, integer *ia, integer *ja, doublereal *a, doublereal *ur, doublereal *vb);
/* Common Block Declarations */
static struct {
integer in;
integer is;
integer isym;
integer itmax;
integer level;
integer nout;
} itcom1_;
#define itcom1_1 itcom1_
static struct {
logical adapt;
logical betadt;
logical caseii;
logical halt;
logical partad;
} itcom2_;
#define itcom2_1 itcom2_
static struct {
doublereal bdelnm;
doublereal betab;
doublereal cme;
doublereal delnnm;
doublereal delsnm;
doublereal ff;
doublereal gamma;
doublereal omega;
doublereal qa;
doublereal qt;
doublereal rho;
doublereal rrr;
doublereal sige;
doublereal sme;
doublereal specr;
doublereal spr;
doublereal drelpr;
doublereal stptst;
doublereal udnm;
doublereal zeta;
} itcom3_;
#define itcom3_1 itcom3_
/* Table of constant values */
static integer c__1 = 1;
static integer c__0 = 0;
static doublereal c_b21 = 0.;
static integer c__2 = 2;
static integer c__3 = 3;
static integer c__4 = 4;
static doublereal c_b286 = 1.;
static integer c__5 = 5;
static integer c__6 = 6;
static integer c__7 = 7;
/* Subroutine */
int jcg_(integer *n, integer *ia, integer *ja, doublereal *a, doublereal *rhs, doublereal *u,
integer *iwksp, integer *nw, doublereal *wksp, integer *iparm, doublereal *rparm, integer *ierr)
{
/* Local variables */
static integer n3, nb, ib1, ib2, ib3, ib4, ib5, ier;
static doublereal tol;
static doublereal temp;
static integer loop;
static doublereal time1, time2;
static real timi1, timj1, timi2, timj2;
static integer idgts;
static doublereal digit1, digit2;
static integer itmax1;
static integer ierper;
/* ITPACK 2C MAIN SUBROUTINE JCG (JACOBI CONJUGATE GRADIENT) */
/* EACH OF THE MAIN SUBROUTINES: */
/* JCG, JSI, SOR, SSORCG, SSORSI, RSCG, RSSI */
/* CAN BE USED INDEPENDENTLY OF THE OTHERS */
/* THIS SUBROUTINE, JCG, DRIVES THE JACOBI CONJUGATE */
/* GRADIENT ALGORITHM. */
/* */
/* ... PARAMETER LIST: */
/* */
/* N INPUT INTEGER. DIMENSION OF THE MATRIX. */
/* IA,JA INPUT INTEGER VECTORS. THE TWO INTEGER ARRAYS OF */
/* THE SPARSE MATRIX REPRESENTATION. */
/* A INPUT D.P. VECTOR. THE D.P. ARRAY OF THE SPARSE */
/* MATRIX REPRESENTATION. */
/* RHS INPUT D.P. VECTOR. CONTAINS THE RIGHT HAND SIDE */
/* OF THE MATRIX PROBLEM. */
/* U INPUT/OUTPUT D.P. VECTOR. ON INPUT, U CONTAINS THE */
/* INITIAL GUESS TO THE SOLUTION. ON OUTPUT, IT CONTAINS */
/* THE LATEST ESTIMATE TO THE SOLUTION. */
/* IWKSP INTEGER VECTOR WORKSPACE OF LENGTH 3*N */
/* NW INPUT INTEGER. LENGTH OF AVAILABLE WKSP. ON OUTPUT, */
/* IPARM(8) IS AMOUNT USED. */
/* WKSP D.P. VECTOR USED FOR WORKING SPACE. JACOBI CONJUGATE */
/* GRADIENT NEEDS THIS TO BE IN LENGTH AT LEAST */
/* 4*N + 2*ITMAX, IF ISYM = 0 (SYMMETRIC STORAGE) */
/* 4*N + 4*ITMAX, IF ISYM = 1 (NONSYMMETRIC STORAGE) */
/* HERE ITMAX = IPARM(1) AND ISYM = IPARM(5) */
/* (ITMAX IS THE MAXIMUM ALLOWABLE NUMBER OF ITERATIONS) */
/* IPARM INTEGER VECTOR OF LENGTH 12. ALLOWS USER TO SPECIFY */
/* SOME INTEGER PARAMETERS WHICH AFFECT THE METHOD. */
/* RPARM D.P. VECTOR OF LENGTH 12. ALLOWS USER TO SPECIFY SOME */
/* D.P. PARAMETERS WHICH AFFECT THE METHOD. */
/* IER OUTPUT INTEGER. ERROR FLAG. (= IERR) */
/* */
/* ... JCG SUBPROGRAM REFERENCES: */
/* */
/* FROM ITPACK BISRCH, CHGCON, DETERM, DFAULT, ECHALL, */
/* ECHOUT, EIGVNS, EIGVSS, EQRT1S, ITERM, TIMER, */
/* ITJCG, IVFILL, PARCON, PERMAT, */
/* PERROR, PERVEC, PJAC, PMULT, PRBNDX, */
/* PSTOP, QSORT, DAXPY, SBELM, SCAL, DCOPY, */
/* DDOT, SUM3, UNSCAL, VEVMW, VFILL, VOUT, */
/* WEVMW, ZBRENT */
/* SYSTEM DABS, DLOG10, DBLE(AMAX0), DMAX1, MOD, DSQRT */
/* */
/* VERSION: ITPACK 2C (MARCH 1982) */
/* */
/* CODE WRITTEN BY: DAVID KINCAID, ROGER GRIMES, JOHN RESPESS */
/* CENTER FOR NUMERICAL ANALYSIS */
/* UNIVERSITY OF TEXAS */
/* AUSTIN, TX 78712 */
/* (512) 471-1242 */
/* */
/* FOR ADDITIONAL DETAILS ON THE */
/* (A) SUBROUTINE SEE TOMS ARTICLE 1982 */
/* (B) ALGORITHM SEE CNA REPORT 150 */
/* */
/* BASED ON THEORY BY: DAVID YOUNG, DAVID KINCAID, LOU HAGEMAN */
/* */
/* REFERENCE THE BOOK: APPLIED ITERATIVE METHODS */
/* L. HAGEMAN, D. YOUNG */
/* ACADEMIC PRESS, 1981 */
/* */
/* ************************************************** */
/* * IMPORTANT NOTE * */
/* * * */
/* * WHEN INSTALLING ITPACK ROUTINES ON A * */
/* * DIFFERENT COMPUTER, RESET SOME OF THE VALUES * */
/* * IN SUBROUTNE DFAULT. MOST IMPORTANT ARE * */
/* * * */
/* * DRELPR MACHINE RELATIVE PRECISION * */
/* * RPARM(1) STOPPING CRITERION * */
/* * * */
/* * ALSO CHANGE SYSTEM-DEPENDENT ROUTINE * */
/* * SECOND USED IN TIMER * */
/* * * */
/* ************************************************** */
/* */
/* ... VARIABLES IN COMMON BLOCK - ITCOM1 */
/* */
/* IN - ITERATION NUMBER */
/* IS - ITERATION NUMBER WHEN PARAMETERS LAST CHANGED */
/* ISYM - SYMMETRIC/NONSYMMETRIC STORAGE FORMAT SWITCH */
/* ITMAX - MAXIMUM NUMBER OF ITERATIONS ALLOWED */
/* LEVEL - LEVEL OF OUTPUT CONTROL SWITCH */
/* NOUT - OUTPUT UNIT NUMBER */
/* */
/* ... VARIABLES IN COMMON BLOCK - ITCOM2 */
/* */
/* ADAPT - FULLY ADAPTIVE PROCEDURE SWITCH */
/* BETADT - SWITCH FOR ADAPTIVE DETERMINATION OF BETA */
/* CASEII - ADAPTIVE PROCEDURE CASE SWITCH */
/* HALT - STOPPING TEST SWITCH */
/* PARTAD - PARTIALLY ADAPTIVE PROCEDURE SWITCH */
/* */
/* ... VARIABLES IN COMMON BLOCK - ITCOM3 */
/* */
/* BDELNM - TWO NORM OF B TIMES DELTA-SUPER-N */
/* BETAB - ESTIMATE FOR THE SPECTRAL RADIUS OF LU MATRIX */
/* CME - ESTIMATE OF LARGEST EIGENVALUE */
/* DELNNM - INNER PRODUCT OF PSEUDO-RESIDUAL AT ITERATION N */
/* DELSNM - INNER PRODUCT OF PSEUDO-RESIDUAL AT ITERATION S */
/* FF - ADAPTIVE PROCEDURE DAMPING FACTOR */
/* GAMMA - ACCELERATION PARAMETER */
/* OMEGA - OVERRELAXATION PARAMETER FOR SOR AND SSOR */
/* QA - PSEUDO-RESIDUAL RATIO */
/* QT - VIRTUAL SPECTRAL RADIUS */
/* RHO - ACCELERATION PARAMETER */
/* RRR - ADAPTIVE PARAMETER */
/* SIGE - PARAMETER SIGMA-SUB-E */
/* SME - ESTIMATE OF SMALLEST EIGENVALUE */
/* SPECR - SPECTRAL RADIUS ESTIMATE FOR SSOR */
/* DRELPR - MACHINE RELATIVE PRECISION */
/* STPTST - STOPPING PARAMETER */
/* UDNM - TWO NORM OF U */
/* ZETA - STOPPING CRITERION */
itcom1_1.level = iparm[1];
itcom1_1.nout = iparm[3];
ier = 0;
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -