0024624: Lost word in license statement in source files
[occt.git] / src / AdvApp2Var / AdvApp2Var_ApproxF2var.cxx
CommitLineData
973c2be1 1// Copyright (c) 1999-2014 OPEN CASCADE SAS
7fd59977 2//
973c2be1 3// This file is part of Open CASCADE Technology software library.
b311480e 4//
d5f74e42 5// This library is free software; you can redistribute it and/or modify it under
6// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 7// by the Free Software Foundation, with special exception defined in the file
8// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
9// distribution for complete text of the license and disclaimer of any warranty.
7fd59977 10//
973c2be1 11// Alternatively, this file may be used under the terms of Open CASCADE
12// commercial license or contractual agreement.
b311480e 13
14// AdvApp2Var_ApproxF2var.cxx
7fd59977 15#include <math.h>
16#include <AdvApp2Var_SysBase.hxx>
17#include <AdvApp2Var_MathBase.hxx>
18#include <AdvApp2Var_Data_f2c.hxx>
19#include <AdvApp2Var_Data.hxx>
20#include <AdvApp2Var_ApproxF2var.hxx>
21
22
23static
24int mmjacpt_(const integer *ndimen,
25 const integer *ncoefu,
26 const integer *ncoefv,
27 const integer *iordru,
28 const integer *iordrv,
29 const doublereal *ptclgd,
30 doublereal *ptcaux,
31 doublereal *ptccan);
32
33
34
35static
36int mma2ce2_(integer *numdec,
37 integer *ndimen,
38 integer *nbsesp,
39 integer *ndimse,
40 integer *ndminu,
41 integer *ndminv,
42 integer *ndguli,
43 integer *ndgvli,
44 integer *ndjacu,
45 integer *ndjacv,
46 integer *iordru,
47 integer *iordrv,
48 integer *nbpntu,
49 integer *nbpntv,
50 doublereal *epsapr,
51 doublereal *sosotb,
52 doublereal *disotb,
53 doublereal *soditb,
54 doublereal *diditb,
55 doublereal *gssutb,
56 doublereal *gssvtb,
57 doublereal *xmaxju,
58 doublereal *xmaxjv,
59 doublereal *vecerr,
60 doublereal *chpair,
61 doublereal *chimpr,
62 doublereal *patjac,
63 doublereal *errmax,
64 doublereal *errmoy,
65 integer *ndegpu,
66 integer *ndegpv,
67 integer *itydec,
68 integer *iercod);
69
70static
71int mma2cfu_(integer *ndujac,
72 integer *nbpntu,
73 integer *nbpntv,
74 doublereal *sosotb,
75 doublereal *disotb,
76 doublereal *soditb,
77 doublereal *diditb,
78 doublereal *gssutb,
79 doublereal *chpair,
80 doublereal *chimpr);
81
82static
83int mma2cfv_(integer *ndvjac,
84 integer *mindgu,
85 integer *maxdgu,
86 integer *nbpntv,
87 doublereal *gssvtb,
88 doublereal *chpair,
89 doublereal *chimpr,
90 doublereal *patjac);
91
92static
93int mma2er1_(integer *ndjacu,
94 integer *ndjacv,
95 integer *ndimen,
96 integer *mindgu,
97 integer *maxdgu,
98 integer *mindgv,
99 integer *maxdgv,
100 integer *iordru,
101 integer *iordrv,
102 doublereal *xmaxju,
103 doublereal *xmaxjv,
104 doublereal *patjac,
105 doublereal *vecerr,
106 doublereal *erreur);
107
108static
109int mma2er2_(integer *ndjacu,
110 integer *ndjacv,
111 integer *ndimen,
112 integer *mindgu,
113 integer *maxdgu,
114 integer *mindgv,
115 integer *maxdgv,
116 integer *iordru,
117 integer *iordrv,
118 doublereal *xmaxju,
119 doublereal *xmaxjv,
120 doublereal *patjac,
121 doublereal *epmscut,
122 doublereal *vecerr,
123 doublereal *erreur,
124 integer *newdgu,
125 integer *newdgv);
126
127static
128int mma2moy_(integer *ndgumx,
129 integer *ndgvmx,
130 integer *ndimen,
131 integer *mindgu,
132 integer *maxdgu,
133 integer *mindgv,
134 integer *maxdgv,
135 integer *iordru,
136 integer *iordrv,
137 doublereal *patjac,
138 doublereal *errmoy);
139
140static
141int mma2ds2_(integer *ndimen,
142 doublereal *uintfn,
143 doublereal *vintfn,
41194117 144 const AdvApp2Var_EvaluatorFunc2Var& foncnp,
7fd59977 145 integer *nbpntu,
146 integer *nbpntv,
147 doublereal *urootb,
148 doublereal *vrootb,
149 integer *iiuouv,
150 doublereal *sosotb,
151 doublereal *disotb,
152 doublereal *soditb,
153 doublereal *diditb,
154 doublereal *fpntab,
155 doublereal *ttable,
156 integer *iercod);
157
158
159
160
161static
162int mma1fdi_(integer *ndimen,
163 doublereal *uvfonc,
41194117 164 const AdvApp2Var_EvaluatorFunc2Var& foncnp,
7fd59977 165 integer *isofav,
166 doublereal *tconst,
167 integer *nbroot,
168 doublereal *ttable,
169 integer *iordre,
170 integer *ideriv,
171 doublereal *fpntab,
172 doublereal *somtab,
173 doublereal *diftab,
174 doublereal *contr1,
175 doublereal *contr2,
176 integer *iercod);
177
178static
179int mma1cdi_(integer *ndimen,
180 integer *nbroot,
181 doublereal *rootlg,
182 integer *iordre,
183 doublereal *contr1,
184 doublereal *contr2,
185 doublereal *somtab,
186 doublereal *diftab,
187 doublereal *fpntab,
188 doublereal *hermit,
189 integer *iercod);
190static
191int mma1jak_(integer *ndimen,
192 integer *nbroot,
193 integer *iordre,
194 integer *ndgjac,
195 doublereal *somtab,
196 doublereal *diftab,
197 doublereal *cgauss,
198 doublereal *crvjac,
199 integer *iercod);
200static
201int mma1cnt_(integer *ndimen,
202 integer *iordre,
203 doublereal *contr1,
204 doublereal *contr2,
205 doublereal *hermit,
206 integer *ndgjac,
207 doublereal *crvjac);
208
209static
210int mma1fer_(integer *ndimen,
211 integer *nbsesp,
212 integer *ndimse,
213 integer *iordre,
214 integer *ndgjac,
215 doublereal *crvjac,
216 integer *ncflim,
217 doublereal *epsapr,
218 doublereal *ycvmax,
219 doublereal *errmax,
220 doublereal *errmoy,
221 integer *ncoeff,
222 integer *iercod);
223
224static
225int mma1noc_(doublereal *dfuvin,
226 integer *ndimen,
227 integer *iordre,
228 doublereal *cntrin,
229 doublereal *duvout,
230 integer *isofav,
231 integer *ideriv,
232 doublereal *cntout);
233
234
235static
236 int mmmapcoe_(integer *ndim,
237 integer *ndgjac,
238 integer *iordre,
239 integer *nbpnts,
240 doublereal *somtab,
241 doublereal *diftab,
242 doublereal *gsstab,
243 doublereal *crvjac);
244
245static
246 int mmaperm_(integer *ncofmx,
247 integer *ndim,
248 integer *ncoeff,
249 integer *iordre,
250 doublereal *crvjac,
251 integer *ncfnew,
252 doublereal *errmoy);
253
254
255#define mmapgss_1 mmapgss_
256#define mmapgs0_1 mmapgs0_
257#define mmapgs1_1 mmapgs1_
258#define mmapgs2_1 mmapgs2_
259
260//=======================================================================
261//function : mma1cdi_
262//purpose :
263//=======================================================================
264int mma1cdi_(integer *ndimen,
265 integer *nbroot,
266 doublereal *rootlg,
267 integer *iordre,
268 doublereal *contr1,
269 doublereal *contr2,
270 doublereal *somtab,
271 doublereal *diftab,
272 doublereal *fpntab,
273 doublereal *hermit,
274 integer *iercod)
275{
1ef32e96 276 integer c__1 = 1;
7fd59977 277
278 /* System generated locals */
279 integer contr1_dim1, contr1_offset, contr2_dim1, contr2_offset,
280 somtab_dim1, somtab_offset, diftab_dim1, diftab_offset,
281 fpntab_dim1, fpntab_offset, hermit_dim1, hermit_offset, i__1,
282 i__2, i__3;
41194117 283
7fd59977 284 /* Local variables */
1ef32e96
RL
285 integer nroo2, ncfhe, nd, ii, kk;
286 integer ibb, kkm, kkp;
1d47d8d0 287 doublereal bid1, bid2, bid3 = 0.;
7fd59977 288
7fd59977 289/* **********************************************************************
290*/
0d969553 291/* FUNCTION : */
7fd59977 292/* ---------- */
0d969553
Y
293/* Discretisation on the parameters of interpolation polynomes */
294/* constraints of order IORDRE. */
7fd59977 295
0d969553 296/* KEYWORDS : */
7fd59977 297/* ----------- */
0d969553 298/* ALL, AB_SPECIFI::CONTRAINTE&, DISCRETISATION, &POINT */
7fd59977 299
0d969553 300/* INPUT ARGUMENTS : */
7fd59977 301/* ------------------ */
0d969553
Y
302/* NDIMEN: Space dimension. */
303/* NBROOT: Number of INTERNAL discretisation parameters. */
304/* It is also the root number Legendre polynome where */
305/* the discretization is performed. */
306/* ROOTLG: Table of discretization parameters ON (-1,1). */
307/* IORDRE: Order of constraint imposed to the extremities of the iso. */
308/* = 0, the extremities of the iso are calculated */
309/* = 1, additionally, the 1st derivative in the direction */
310/* of the iso is calculated. */
311/* = 2, additionally, the 2nd derivative in the direction */
312/* of the iso is calculated. */
313/* CONTR1: Contains, if IORDRE>=0, values IORDRE+1 in TTABLE(0)
314*/
315/* (1st extremity) of derivatives of F(Uc,Ve) or F(Ue,Vc), */
316/* see below. */
317/* CONTR2: Contains, if IORDRE>=0, values IORDRE+1 in */
318/* TTABLE(NBROOT+1) (2nd extremity) of: */
319/* If ISOFAV=1, derived of order IDERIV by U, derived */
320/* ordre 0 to IORDRE by V of F(Uc,Ve) or Uc=TCONST */
321/* (fixed iso value) and Ve is the fixed extremity. */
322/* If ISOFAV=2, derivative of order IDERIV by V, derivative */
323/* of order 0 to IORDRE by U of F(Ue,Vc) or Vc=TCONST */
324/* (fixed iso value) and Ue is the fixed extremity. */
325
326/* SOMTAB: Table of NBROOT/2 sums of 2 index points */
327/* NBROOT-II+1 and II, for II = 1, NBROOT/2. */
328/* DIFTAB: Table of NBROOT/2 differences of 2 index points */
329/* NBROOT-II+1 and II, for II = 1, NBROOT/2. */
330
331/* OUTPUT ARGUMENTS : */
7fd59977 332/* ------------------- */
0d969553
Y
333/* SOMTAB: Table of NBROOT/2 sums of 2 index points */
334/* NBROOT-II+1 and II, for II = 1, NBROOT/2 */
335/* DIFTAB: Table of NBROOT/2 differences of 2 index points */
336/* NBROOT-II+1 and II, for II = 1, NBROOT/2 */
337/* FPNTAB: Auxiliary table. */
338/* HERMIT: Table of coeff. 2*(IORDRE+1) Hermite polynoms */
339/* of degree 2*IORDRE+1. */
340/* IERCOD: Error code, */
341/* = 0, Everythig is OK */
342/* = 1, The value of IORDRE is out of (0,2) */
343/* COMMON USED : */
7fd59977 344/* ---------------- */
345
0d969553 346/* REFERENCES CALLED : */
7fd59977 347/* ----------------------- */
348
0d969553 349/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 350/* ----------------------------------- */
0d969553
Y
351/* The results of discretization are arranged in 2 tables */
352/* SOMTAB and DIFTAB to earn time during the */
353/* calculation of coefficients of the approximation curve. */
7fd59977 354
0d969553
Y
355/* If NBROOT is uneven in SOMTAB(0,*) and DIFTAB(0,*) one stores */
356/* the values of the median root of Legendre (0.D0 in (-1,1)). */
7fd59977 357
7fd59977 358/* **********************************************************************
359*/
360
0d969553 361/* Name of the routine */
7fd59977 362
363
364 /* Parameter adjustments */
365 diftab_dim1 = *nbroot / 2 + 1;
366 diftab_offset = diftab_dim1;
367 diftab -= diftab_offset;
368 somtab_dim1 = *nbroot / 2 + 1;
369 somtab_offset = somtab_dim1;
370 somtab -= somtab_offset;
371 --rootlg;
372 hermit_dim1 = (*iordre << 1) + 2;
373 hermit_offset = hermit_dim1;
374 hermit -= hermit_offset;
375 fpntab_dim1 = *nbroot;
376 fpntab_offset = fpntab_dim1 + 1;
377 fpntab -= fpntab_offset;
378 contr2_dim1 = *ndimen;
379 contr2_offset = contr2_dim1 + 1;
380 contr2 -= contr2_offset;
381 contr1_dim1 = *ndimen;
382 contr1_offset = contr1_dim1 + 1;
383 contr1 -= contr1_offset;
384
385 /* Function Body */
386 ibb = AdvApp2Var_SysBase::mnfndeb_();
387 if (ibb >= 3) {
388 AdvApp2Var_SysBase::mgenmsg_("MMA1CDI", 7L);
389 }
390 *iercod = 0;
391
0d969553 392/* --- Recuperate 2*(IORDRE+1) coeff of 2*(IORDRE+1) of Hermite polynom ---
7fd59977 393*/
394
395 AdvApp2Var_ApproxF2var::mma1her_(iordre, &hermit[hermit_offset], iercod);
396 if (*iercod > 0) {
397 goto L9100;
398 }
399
0d969553 400/* ------------------- Discretization of Hermite polynoms -----------
7fd59977 401*/
402
403 ncfhe = (*iordre + 1) << 1;
404 i__1 = ncfhe;
405 for (ii = 1; ii <= i__1; ++ii) {
406 i__2 = *nbroot;
407 for (kk = 1; kk <= i__2; ++kk) {
408 AdvApp2Var_MathBase::mmmpocur_(&ncfhe, &c__1, &ncfhe, &hermit[ii * hermit_dim1], &
409 rootlg[kk], &fpntab[kk + ii * fpntab_dim1]);
410/* L200: */
411 }
412/* L100: */
413 }
414
0d969553 415/* ---- Discretizations of boundary polynoms are taken ----
7fd59977 416*/
417
418 nroo2 = *nbroot / 2;
419 i__1 = *ndimen;
420 for (nd = 1; nd <= i__1; ++nd) {
421 i__2 = *iordre + 1;
422 for (ii = 1; ii <= i__2; ++ii) {
423 bid1 = contr1[nd + ii * contr1_dim1];
424 bid2 = contr2[nd + ii * contr2_dim1];
425 i__3 = nroo2;
426 for (kk = 1; kk <= i__3; ++kk) {
427 kkm = nroo2 - kk + 1;
428 bid3 = bid1 * fpntab[kkm + ((ii << 1) - 1) * fpntab_dim1] +
429 bid2 * fpntab[kkm + (ii << 1) * fpntab_dim1];
430 somtab[kk + nd * somtab_dim1] -= bid3;
431 diftab[kk + nd * diftab_dim1] += bid3;
432/* L500: */
433 }
434 i__3 = nroo2;
435 for (kk = 1; kk <= i__3; ++kk) {
436 kkp = (*nbroot + 1) / 2 + kk;
437 bid3 = bid1 * fpntab[kkp + ((ii << 1) - 1) * fpntab_dim1] +
438 bid2 * fpntab[kkp + (ii << 1) * fpntab_dim1];
439 somtab[kk + nd * somtab_dim1] -= bid3;
440 diftab[kk + nd * diftab_dim1] -= bid3;
441/* L600: */
442 }
443/* L400: */
444 }
445/* L300: */
446 }
447
0d969553 448/* ------------ Cas when discretization is done on the roots of a -----------
7fd59977 449*/
0d969553 450/* ---------- Legendre polynom of uneven degree, 0 is root --------
7fd59977 451*/
452
453 if (*nbroot % 2 == 1) {
454 i__1 = *ndimen;
455 for (nd = 1; nd <= i__1; ++nd) {
456 i__2 = *iordre + 1;
457 for (ii = 1; ii <= i__2; ++ii) {
458 bid3 = fpntab[nroo2 + 1 + ((ii << 1) - 1) * fpntab_dim1] *
459 contr1[nd + ii * contr1_dim1] + fpntab[nroo2 + 1 + (
460 ii << 1) * fpntab_dim1] * contr2[nd + ii *
461 contr2_dim1];
462/* L800: */
463 }
464 somtab[nd * somtab_dim1] -= bid3;
465 diftab[nd * diftab_dim1] -= bid3;
466/* L700: */
467 }
468 }
469
470 goto L9999;
471
472/* ------------------------------ The End -------------------------------
473*/
0d969553 474/* --> IORDRE is not in the authorized zone. */
7fd59977 475L9100:
476 *iercod = 1;
477 goto L9999;
478
479L9999:
480 if (ibb >= 3) {
481 AdvApp2Var_SysBase::mgsomsg_("MMA1CDI", 7L);
482 }
483 return 0;
484} /* mma1cdi_ */
485
486//=======================================================================
487//function : mma1cnt_
488//purpose :
489//=======================================================================
490int mma1cnt_(integer *ndimen,
491 integer *iordre,
492 doublereal *contr1,
493 doublereal *contr2,
494 doublereal *hermit,
495 integer *ndgjac,
496 doublereal *crvjac)
497{
498 /* System generated locals */
499 integer contr1_dim1, contr1_offset, contr2_dim1, contr2_offset,
500 hermit_dim1, hermit_offset, crvjac_dim1, crvjac_offset, i__1,
501 i__2, i__3;
41194117 502
7fd59977 503 /* Local variables */
1ef32e96
RL
504 integer nd, ii, jj, ibb;
505 doublereal bid;
41194117
K
506
507
7fd59977 508 /* ***********************************************************************
509 */
510
0d969553 511 /* FUNCTION : */
7fd59977 512 /* ---------- */
0d969553 513 /* Add constraint to polynom. */
7fd59977 514
515 /* MOTS CLES : */
516 /* ----------- */
0d969553 517 /* ALL,AB_SPECIFI::COURE&,APPROXIMATION,ADDITION,&CONSTRAINT */
7fd59977 518
0d969553 519 /* INPUT ARGUMENTS : */
7fd59977 520 /* -------------------- */
0d969553
Y
521 /* NDIMEN: Dimension of the space */
522 /* IORDRE: Order of constraint. */
523 /* CONTR1: pt of constraint in -1, from order 0 to IORDRE. */
524 /* CONTR2: Pt of constraint in +1, from order 0 to IORDRE. */
525 /* HERMIT: Table of Hermit polynoms of order IORDRE. */
526 /* CRVJAV: Curve of approximation in Jacobi base. */
7fd59977 527
0d969553 528 /* OUTPUT ARGUMENTS : */
7fd59977 529 /* --------------------- */
0d969553
Y
530 /* CRVJAV: Curve of approximation in Jacobi base */
531 /* to which the polynom of interpolation of constraints is added. */
7fd59977 532
0d969553 533 /* COMMON USED : */
7fd59977 534 /* ------------------ */
535
536
0d969553 537 /* REFERENCES CALLED : */
7fd59977 538 /* --------------------- */
539
540
0d969553 541/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 542/* ----------------------------------- */
543
7fd59977 544/* > */
545/* ***********************************************************************
546 */
547/* DECLARATIONS */
548/* ***********************************************************************
549 */
0d969553 550/* Name of the routine */
7fd59977 551
552/* ***********************************************************************
553 */
554/* INITIALISATIONS */
555/* ***********************************************************************
556 */
557
558 /* Parameter adjustments */
559 hermit_dim1 = (*iordre << 1) + 2;
560 hermit_offset = hermit_dim1;
561 hermit -= hermit_offset;
562 contr2_dim1 = *ndimen;
563 contr2_offset = contr2_dim1 + 1;
564 contr2 -= contr2_offset;
565 contr1_dim1 = *ndimen;
566 contr1_offset = contr1_dim1 + 1;
567 contr1 -= contr1_offset;
568 crvjac_dim1 = *ndgjac + 1;
569 crvjac_offset = crvjac_dim1;
570 crvjac -= crvjac_offset;
571
572 /* Function Body */
573 ibb = AdvApp2Var_SysBase::mnfndeb_();
574 if (ibb >= 3) {
575 AdvApp2Var_SysBase::mgenmsg_("MMA1CNT", 7L);
576 }
577
578/* ***********************************************************************
579 */
0d969553 580/* Processing */
7fd59977 581/* ***********************************************************************
582 */
583
584 i__1 = *ndimen;
585 for (nd = 1; nd <= i__1; ++nd) {
586 i__2 = (*iordre << 1) + 1;
587 for (ii = 0; ii <= i__2; ++ii) {
588 bid = 0.;
589 i__3 = *iordre + 1;
590 for (jj = 1; jj <= i__3; ++jj) {
591 bid = bid + contr1[nd + jj * contr1_dim1] *
592 hermit[ii + ((jj << 1) - 1) * hermit_dim1] +
593 contr2[nd + jj * contr2_dim1] * hermit[ii + (jj << 1) * hermit_dim1];
594 /* L300: */
595 }
596 crvjac[ii + nd * crvjac_dim1] = bid;
597 /* L200: */
598 }
599 /* L100: */
600 }
601
602/* ***********************************************************************
603 */
0d969553 604/* RETURN CALLING PROGRAM */
7fd59977 605/* ***********************************************************************
606 */
607
608 if (ibb >= 3) {
609 AdvApp2Var_SysBase::mgsomsg_("MMA1CNT", 7L);
610 }
611
612 return 0 ;
613} /* mma1cnt_ */
614
615//=======================================================================
616//function : mma1fdi_
617//purpose :
618//=======================================================================
619int mma1fdi_(integer *ndimen,
620 doublereal *uvfonc,
41194117 621 const AdvApp2Var_EvaluatorFunc2Var& foncnp,
7fd59977 622 integer *isofav,
623 doublereal *tconst,
624 integer *nbroot,
625 doublereal *ttable,
626 integer *iordre,
627 integer *ideriv,
628 doublereal *fpntab,
629 doublereal *somtab,
630 doublereal *diftab,
631 doublereal *contr1,
632 doublereal *contr2,
633 integer *iercod)
634{
635 /* System generated locals */
636 integer fpntab_dim1, somtab_dim1, somtab_offset, diftab_dim1,
637 diftab_offset, contr1_dim1, contr1_offset, contr2_dim1,
638 contr2_offset, i__1, i__2;
639 doublereal d__1;
41194117 640
7fd59977 641 /* Local variables */
1ef32e96
RL
642 integer ideb, ifin, nroo2, ideru, iderv;
643 doublereal renor;
644 integer ii, nd, ibb, iim, nbp, iip;
645 doublereal bid1, bid2;
41194117 646
7fd59977 647/* **********************************************************************
648*/
649
0d969553 650/* FUNCTION : */
7fd59977 651/* ---------- */
0d969553
Y
652/* DiscretiZation of a non-polynomial function F(U,V) or of */
653/* its derivative with fixed isoparameter. */
7fd59977 654
0d969553 655/* KEYWORDS : */
7fd59977 656/* ----------- */
0d969553 657/* ALL, AB_SPECIFI::FONCTION&, DISCRETISATION, &POINT */
7fd59977 658
0d969553 659/* INPUT ARGUMENTS : */
7fd59977 660/* ------------------ */
0d969553
Y
661/* NDIMEN: Space dimension. */
662/* UVFONC: Limits of the path of definition by U and by V of the approximated function */
663/* FONCNP: The NAME of the non-polynomial function to be approximated */
664/* (external program). */
665/* ISOFAV: Fixed isoparameter for the discretization; */
666/* = 1, discretization with fixed U and variable V. */
667/* = 2, discretization with fixed V and variable U. */
668/* TCONST: Iso value is also fixed. */
669/* NBROOT: Number of INTERNAL discretization parameters. */
670/* (if there are constraints, 2 extremities should be added).
671*/
672/* This is also the root number of the Legendre polynom where */
673/* the discretization is done. */
674/* TTABLE: Table of discretization parameters and of 2 extremities */
675/* (Respectively (-1, NBROOT Legendre roots,1) */
676/* reframed within the adequate interval. */
677/* IORDRE: Order of constraint imposed on the extremities of the iso. */
678/* (If Iso-U, it is necessary to calculate the derivatives by V and vice */
7fd59977 679/* versa). */
0d969553
Y
680/* = 0, the extremities of the iso are calculated. */
681/* = 1, additionally the 1st derivative in the direction of the iso is calculated */
682/* = 2, additionally the 2nd derivative in the direction of the iso is calculated */
683/* IDERIV: Order of derivative transversal to fixed iso (If Iso-U=Uc */
684/* is fixed, the derivative of order IDERIV is discretized by U of */
685/* F(Uc,v). Same if iso-V is fixed). */
686/* Varies from 0 (positioning) to 2 (2nd derivative). */
687
688/* OUTPUT ARGUMENTS : */
7fd59977 689/* ------------------- */
0d969553
Y
690/* FPNTAB: Auxiliary table.
691 SOMTAB: Table of NBROOT/2 sums of 2 index points */
692/* NBROOT-II+1 and II, for II = 1, NBROOT/2 */
693/* DIFTAB: Table of NBROOT/2 differences of 2 index points */
694/* NBROOT-II+1 and II, for II = 1, NBROOT/2 */
695/* CONTR1: Contains, if IORDRE>=0, values IORDRE+1 in TTABLE(0)
696*/
697/* (1st extremity) of derivatives of F(Uc,Ve) or F(Ue,Vc), */
698/* see below. */
699/* CONTR2: Contains, if IORDRE>=0, values IORDRE+1 in */
700/* TTABLE(NBROOT+1) (2nd extremity) of: */
701/* If ISOFAV=1, derived of order IDERIV by U, derived */
702/* ordre 0 to IORDRE by V of F(Uc,Ve) or Uc=TCONST */
703/* (fixed iso value) and Ve is the fixed extremity. */
704/* If ISOFAV=2, derivative of order IDERIV by V, derivative */
705/* of order 0 to IORDRE by U of F(Ue,Vc) or Vc=TCONST */
706/* (fixed iso value) and Ue is the fixed extremity. */
707/* IERCOD: Error code > 100; Pb in evaluation of FONCNP, */
708/* the returned error code is equal to error code of FONCNP + 100. */
709
710/* COMMONS USED : */
7fd59977 711/* ---------------- */
712
0d969553 713/* REFERENCES CALLED : */
7fd59977 714/* ----------------------- */
715
0d969553 716/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 717/* ----------------------------------- */
0d969553
Y
718/* The results of discretization are arranged in 2 tables */
719/* SOMTAB and DIFTAB to earn time during the */
720/* calculation of coefficients of the approximation curve. */
7fd59977 721
0d969553
Y
722/* If NBROOT is uneven in SOMTAB(0,*) and DIFTAB(0,*) one stores */
723/* the values of the median root of Legendre (0.D0 in (-1,1)). */
7fd59977 724
0d969553
Y
725/* Function F(u,v) defined in UVFONC is reparameterized in */
726/* (-1,1)x(-1,1). Then 1st and 2nd derivatives are renormalized. */
7fd59977 727
7fd59977 728/* > */
729/* **********************************************************************
730*/
731
0d969553 732/* Name of the routine */
7fd59977 733
734
735 /* Parameter adjustments */
736 uvfonc -= 3;
737 diftab_dim1 = *nbroot / 2 + 1;
738 diftab_offset = diftab_dim1;
739 diftab -= diftab_offset;
740 somtab_dim1 = *nbroot / 2 + 1;
741 somtab_offset = somtab_dim1;
742 somtab -= somtab_offset;
743 fpntab_dim1 = *ndimen;
744 --fpntab;
745 contr2_dim1 = *ndimen;
746 contr2_offset = contr2_dim1 + 1;
747 contr2 -= contr2_offset;
748 contr1_dim1 = *ndimen;
749 contr1_offset = contr1_dim1 + 1;
750 contr1 -= contr1_offset;
751
752 /* Function Body */
753 ibb = AdvApp2Var_SysBase::mnfndeb_();
754 if (ibb >= 3) {
755 AdvApp2Var_SysBase::mgenmsg_("MMA1FDI", 7L);
756 }
757 *iercod = 0;
758
0d969553 759/* --------------- Definition of the nb of points to calculate --------------
7fd59977 760*/
0d969553 761/* --> If constraints, the limits are also taken */
7fd59977 762 if (*iordre >= 0) {
763 ideb = 0;
764 ifin = *nbroot + 1;
0d969553 765/* --> Otherwise, only Legendre roots (reframed) are used
7fd59977 766. */
767 } else {
768 ideb = 1;
769 ifin = *nbroot;
770 }
0d969553 771/* --> Nb of point to calculate. */
7fd59977 772 nbp = ifin - ideb + 1;
773 nroo2 = *nbroot / 2;
774
0d969553 775/* --------------- Determination of the order of global derivation --------
7fd59977 776*/
0d969553
Y
777/* --> ISOFAV takes only values 1 or 2. */
778/* if Iso-U, derive by U of order IDERIV */
7fd59977 779 if (*isofav == 1) {
780 ideru = *ideriv;
781 iderv = 0;
782 d__1 = (uvfonc[4] - uvfonc[3]) / 2.;
783 renor = AdvApp2Var_MathBase::pow__di(&d__1, ideriv);
0d969553 784/* if Iso-V, derive by V of order IDERIV */
7fd59977 785 } else {
786 ideru = 0;
787 iderv = *ideriv;
788 d__1 = (uvfonc[6] - uvfonc[5]) / 2.;
789 renor = AdvApp2Var_MathBase::pow__di(&d__1, ideriv);
790 }
791
0d969553 792/* ----------- Discretization on roots of the ---------------
7fd59977 793*/
0d969553 794/* ---------------------- Legendre polynom of degree NBROOT -------------------
7fd59977 795*/
796
fadcea2c 797 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (ndimen,
7fd59977 798 &uvfonc[3],
799 &uvfonc[5],
800 isofav,
801 tconst,
802 &nbp,
803 &ttable[ideb],
804 &ideru,
805 &iderv,
806 &fpntab[ideb * fpntab_dim1 + 1],
807 iercod);
808 if (*iercod > 0) {
809 goto L9999;
810 }
811 i__1 = *ndimen;
812 for (nd = 1; nd <= i__1; ++nd) {
813 i__2 = nroo2;
814 for (ii = 1; ii <= i__2; ++ii) {
815 iip = (*nbroot + 1) / 2 + ii;
816 iim = nroo2 - ii + 1;
817 bid1 = fpntab[nd + iim * fpntab_dim1];
818 bid2 = fpntab[nd + iip * fpntab_dim1];
819 somtab[ii + nd * somtab_dim1] = renor * (bid2 + bid1);
820 diftab[ii + nd * diftab_dim1] = renor * (bid2 - bid1);
821/* L200: */
822 }
823/* L100: */
824 }
825
0d969553 826/* ------------ Case when discretisation is done on roots of a ----
7fd59977 827*/
0d969553 828/* ---------- Legendre polynom of uneven degree, 0 is root --------
7fd59977 829*/
830
831 if (*nbroot % 2 == 1) {
832 i__1 = *ndimen;
833 for (nd = 1; nd <= i__1; ++nd) {
834 somtab[nd * somtab_dim1] = renor * fpntab[nd + (nroo2 + 1) *
835 fpntab_dim1];
836 diftab[nd * diftab_dim1] = renor * fpntab[nd + (nroo2 + 1) *
837 fpntab_dim1];
838/* L300: */
839 }
840 } else {
841 i__1 = *ndimen;
842 for (nd = 1; nd <= i__1; ++nd) {
843 somtab[nd * somtab_dim1] = 0.;
844 diftab[nd * diftab_dim1] = 0.;
845 }
846 }
847
848
0d969553 849/* --------------------- Take into account constraints ----------------
7fd59977 850*/
851
852 if (*iordre >= 0) {
0d969553 853/* --> Recover already calculated extremities. */
7fd59977 854 i__1 = *ndimen;
855 for (nd = 1; nd <= i__1; ++nd) {
856 contr1[nd + contr1_dim1] = renor * fpntab[nd];
857 contr2[nd + contr2_dim1] = renor * fpntab[nd + (*nbroot + 1) *
858 fpntab_dim1];
859/* L400: */
860 }
0d969553 861/* --> Nb of points to calculate/call to FONCNP */
7fd59977 862 nbp = 1;
0d969553 863/* If Iso-U, derive by V till order IORDRE */
7fd59977 864 if (*isofav == 1) {
0d969553 865/* --> Factor of normalisation 1st derivative. */
7fd59977 866 bid1 = (uvfonc[6] - uvfonc[5]) / 2.;
867 i__1 = *iordre;
868 for (iderv = 1; iderv <= i__1; ++iderv) {
fadcea2c
RL
869 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (
870 ndimen, &uvfonc[3], &uvfonc[5], isofav, tconst, &
871 nbp, ttable, &ideru, &iderv, &contr1[(iderv + 1) *
7fd59977 872 contr1_dim1 + 1], iercod);
873 if (*iercod > 0) {
874 goto L9999;
875 }
876/* L500: */
877 }
878 i__1 = *iordre;
879 for (iderv = 1; iderv <= i__1; ++iderv) {
fadcea2c
RL
880 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (
881 ndimen, &uvfonc[3], &uvfonc[5], isofav, tconst, &
882 nbp, &ttable[*nbroot + 1], &ideru, &iderv, &contr2[(
7fd59977 883 iderv + 1) * contr2_dim1 + 1], iercod);
884 if (*iercod > 0) {
885 goto L9999;
886 }
887/* L510: */
888 }
0d969553 889/* If Iso-V, derive by U till order IORDRE */
7fd59977 890 } else {
0d969553 891/* --> Factor of normalization 1st derivative. */
7fd59977 892 bid1 = (uvfonc[4] - uvfonc[3]) / 2.;
893 i__1 = *iordre;
894 for (ideru = 1; ideru <= i__1; ++ideru) {
fadcea2c
RL
895 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (
896 ndimen, &uvfonc[3], &uvfonc[5], isofav, tconst, &
897 nbp, ttable, &ideru, &iderv, &contr1[(ideru + 1) *
7fd59977 898 contr1_dim1 + 1], iercod);
899 if (*iercod > 0) {
900 goto L9999;
901 }
902/* L600: */
903 }
904 i__1 = *iordre;
905 for (ideru = 1; ideru <= i__1; ++ideru) {
fadcea2c
RL
906 (*const_cast <AdvApp2Var_EvaluatorFunc2Var*> (&foncnp)).Evaluate (
907 ndimen, &uvfonc[3], &uvfonc[5], isofav, tconst, &
908 nbp, &ttable[*nbroot + 1], &ideru, &iderv, &contr2[(
7fd59977 909 ideru + 1) * contr2_dim1 + 1], iercod);
910 if (*iercod > 0) {
911 goto L9999;
912 }
913/* L610: */
914 }
915 }
916
0d969553 917/* ------------------------- Normalization of derivatives -------------
7fd59977 918---- */
0d969553 919/* (The function is redefined on (-1,1)*(-1,1)) */
7fd59977 920 bid2 = renor;
921 i__1 = *iordre;
922 for (ii = 1; ii <= i__1; ++ii) {
923 bid2 = bid1 * bid2;
924 i__2 = *ndimen;
925 for (nd = 1; nd <= i__2; ++nd) {
926 contr1[nd + (ii + 1) * contr1_dim1] *= bid2;
927 contr2[nd + (ii + 1) * contr2_dim1] *= bid2;
928/* L710: */
929 }
930/* L700: */
931 }
932 }
933
934/* ------------------------------ The end -------------------------------
935*/
936
937L9999:
938 if (*iercod > 0) {
939 *iercod += 100;
940 AdvApp2Var_SysBase::maermsg_("MMA1FDI", iercod, 7L);
941 }
942 if (ibb >= 3) {
943 AdvApp2Var_SysBase::mgsomsg_("MMA1FDI", 7L);
944 }
945 return 0;
946} /* mma1fdi_ */
947
948//=======================================================================
949//function : mma1fer_
950//purpose :
951//=======================================================================
952int mma1fer_(integer *,//ndimen,
953 integer *nbsesp,
954 integer *ndimse,
955 integer *iordre,
956 integer *ndgjac,
957 doublereal *crvjac,
958 integer *ncflim,
959 doublereal *epsapr,
960 doublereal *ycvmax,
961 doublereal *errmax,
962 doublereal *errmoy,
963 integer *ncoeff,
964 integer *iercod)
965{
966 /* System generated locals */
967 integer crvjac_dim1, crvjac_offset, i__1, i__2;
41194117 968
7fd59977 969 /* Local variables */
1d47d8d0 970 integer idim, ncfja, ncfnw, ndses, ii, kk, ibb, ier;
1ef32e96 971 integer nbr0;
41194117
K
972
973
7fd59977 974/* ***********************************************************************
975 */
976
0d969553 977/* FUNCTION : */
7fd59977 978/* ---------- */
0d969553 979/* Calculate the degree and the errors of approximation of a border. */
7fd59977 980
0d969553 981/* KEYWORDS : */
7fd59977 982/* ----------- */
983/* TOUS,AB_SPECIFI :: COURBE&,TRONCATURE, &PRECISION */
984
0d969553 985/* INPUT ARGUMENTS : */
7fd59977 986/* -------------------- */
7fd59977 987
0d969553
Y
988/* NDIMEN: Total Dimension of the space (sum of dimensions of sub-spaces) */
989/* NBSESP: Number of "independent" sub-spaces. */
990/* NDIMSE: Table of dimensions of sub-spaces. */
991/* IORDRE: Order of constraint at the extremities of the border */
992/* -1 = no constraints, */
993/* 0 = constraints of passage to limits (i.e. C0), */
994/* 1 = C0 + constraintes of 1st derivatives (i.e. C1), */
995/* 2 = C1 + constraintes of 2nd derivatives (i.e. C2). */
258ff83b 996/* NDGJAC: Degree of development in series to use for the calculation */
0d969553
Y
997/* in the base of Jacobi. */
998/* CRVJAC: Table of coeff. of the curve of approximation in the */
999/* base of Jacobi. */
1000/* NCFLIM: Max number of coeff of the polynomial curve */
1001/* of approximation (should be above or equal to */
1002/* 2*IORDRE+2 and below or equal to 50). */
1003/* EPSAPR: Table of errors of approximations that cannot be passed, */
1004/* sub-space by sub-space. */
1005
1006/* OUTPUT ARGUMENTS : */
7fd59977 1007/* --------------------- */
0d969553
Y
1008/* YCVMAX: Auxiliary Table. */
1009/* ERRMAX: Table of errors (sub-space by sub-space) */
1010/* MAXIMUM made in the approximation of FONCNP by */
7fd59977 1011/* COURBE. */
0d969553
Y
1012/* ERRMOY: Table of errors (sub-space by sub-space) */
1013/* AVERAGE made in the approximation of FONCNP by */
7fd59977 1014/* COURBE. */
0d969553
Y
1015/* NCOEFF: Number of significative coeffs. of the calculated "curve". */
1016/* IERCOD: Error code */
7fd59977 1017/* = 0, ok, */
0d969553
Y
1018/* =-1, warning, required tolerance can't be */
1019/* met with coefficients NFCLIM. */
1020/* = 1, order of constraints (IORDRE) is not within authorised values */
258ff83b 1021
7fd59977 1022
0d969553 1023/* COMMONS USED : */
7fd59977 1024/* ------------------ */
1025
0d969553 1026/* REFERENCES CALLED : */
7fd59977 1027/* --------------------- */
1028
0d969553 1029/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 1030/* ----------------------------------- */
7fd59977 1031/* > */
1032/* **********************************************************************
1033*/
1034
0d969553 1035/* Name of the routine */
7fd59977 1036
1037
1038 /* Parameter adjustments */
1039 --ycvmax;
1040 --errmoy;
1041 --errmax;
1042 --epsapr;
1043 --ndimse;
1044 crvjac_dim1 = *ndgjac + 1;
1045 crvjac_offset = crvjac_dim1;
1046 crvjac -= crvjac_offset;
1047
1048 /* Function Body */
1049 ibb = AdvApp2Var_SysBase::mnfndeb_();
1050 if (ibb >= 3) {
1051 AdvApp2Var_SysBase::mgenmsg_("MMA1FER", 7L);
1052 }
1053 *iercod = 0;
1054 idim = 1;
1055 *ncoeff = 0;
1056 ncfja = *ndgjac + 1;
1057
0d969553 1058/* ------------ Calculate the degree of the curve and of the Max error --------
7fd59977 1059*/
0d969553 1060/* -------------- of approximation for all sub-spaces --------
7fd59977 1061*/
1062
1063 i__1 = *nbsesp;
1064 for (ii = 1; ii <= i__1; ++ii) {
1065 ndses = ndimse[ii];
1066
0d969553 1067/* ------------ cutting of coeff. and calculation of Max error -------
7fd59977 1068---- */
1069
1070 AdvApp2Var_MathBase::mmtrpjj_(&ncfja, &ndses, &ncfja, &epsapr[ii], iordre, &crvjac[idim *
1071 crvjac_dim1], &ycvmax[1], &errmax[ii], &ncfnw);
1072
1073/* ******************************************************************
1074**** */
0d969553 1075/* ------------- If precision OK, calculate the average error -------
7fd59977 1076---- */
1077/* ******************************************************************
1078**** */
1079
1080 if (ncfnw <= *ncflim) {
1081 mmaperm_(&ncfja, &ndses, &ncfja, iordre, &crvjac[idim *
1082 crvjac_dim1], &ncfnw, &errmoy[ii]);
41194117 1083 *ncoeff = advapp_max(ncfnw,*ncoeff);
7fd59977 1084
0d969553 1085/* ------------- Set the declined coefficients to 0.D0 -----------
7fd59977 1086-------- */
1087
1088 nbr0 = *ncflim - ncfnw;
1089 if (nbr0 > 0) {
1090 i__2 = ndses;
1091 for (kk = 1; kk <= i__2; ++kk) {
1092 AdvApp2Var_SysBase::mvriraz_(&nbr0,
fadcea2c 1093 &crvjac[ncfnw + (idim + kk - 1) * crvjac_dim1]);
7fd59977 1094/* L200: */
1095 }
1096 }
1097 } else {
1098
1099/* **************************************************************
1100******** */
0d969553 1101/* ------------------- If required precision can't be reached----
7fd59977 1102-------- */
1103/* **************************************************************
1104******** */
1105
1106 *iercod = -1;
1107
0d969553 1108/* ------------------------- calculate the Max error ------------
7fd59977 1109-------- */
1110
1111 AdvApp2Var_MathBase::mmaperx_(&ncfja, &ndses, &ncfja, iordre, &crvjac[idim *
1112 crvjac_dim1], ncflim, &ycvmax[1], &errmax[ii], &ier);
1113 if (ier > 0) {
1114 goto L9100;
1115 }
1116
0d969553 1117/* -------------------- nb of coeff to be returned -------------
7fd59977 1118-------- */
1119
1120 *ncoeff = *ncflim;
1121
0d969553 1122/* ------------------- and calculation of the average error ----
7fd59977 1123-------- */
1124
1125 mmaperm_(&ncfja, &ndses, &ncfja, iordre, &crvjac[idim *
1126 crvjac_dim1], ncflim, &errmoy[ii]);
1127 }
1128 idim += ndses;
1129/* L100: */
1130 }
1131
1132 goto L9999;
1133
1134/* ------------------------------ The end -------------------------------
1135*/
0d969553 1136/* --> The order of constraints is not within autorized values. */
7fd59977 1137L9100:
1138 *iercod = 1;
1139 goto L9999;
1140
1141L9999:
1142 if (*iercod != 0) {
1143 AdvApp2Var_SysBase::maermsg_("MMA1FER", iercod, 7L);
1144 }
1145 if (ibb >= 3) {
1146 AdvApp2Var_SysBase::mgsomsg_("MMA1FER", 7L);
1147 }
1148 return 0;
1149} /* mma1fer_ */
1150
1151
1152//=======================================================================
1153//function : mma1her_
1154//purpose :
1155//=======================================================================
1156int AdvApp2Var_ApproxF2var::mma1her_(const integer *iordre,
1157 doublereal *hermit,
1158 integer *iercod)
1159{
1160 /* System generated locals */
1161 integer hermit_dim1, hermit_offset;
41194117 1162
7fd59977 1163 /* Local variables */
1ef32e96 1164 integer ibb;
41194117 1165
7fd59977 1166
1167
1168/* **********************************************************************
1169*/
1170
0d969553 1171/* FUNCTION : */
7fd59977 1172/* ---------- */
0d969553
Y
1173/* Calculate 2*(IORDRE+1) Hermit polynoms of degree 2*IORDRE+1 */
1174/* on (-1,1) */
7fd59977 1175
0d969553 1176/* KEYWORDS : */
7fd59977 1177/* ----------- */
0d969553 1178/* ALL, AB_SPECIFI::CONTRAINTE&, INTERPOLATION, &POLYNOME */
7fd59977 1179
0d969553 1180/* INPUT ARGUMENTS : */
7fd59977 1181/* ------------------ */
0d969553
Y
1182/* IORDRE: Order of constraint. */
1183/* = 0, Polynom of interpolation of order C0 on (-1,1). */
1184/* = 1, Polynom of interpolation of order C0 and C1 on (-1,1). */
1185/* = 2, Polynom of interpolation of order C0, C1 and C2 on (-1,1).
7fd59977 1186*/
1187
0d969553 1188/* OUTPUT ARGUMENTS : */
7fd59977 1189/* ------------------- */
0d969553
Y
1190/* HERMIT: Table of 2*IORDRE+2 coeff. of each of 2*(IORDRE+1) */
1191/* HERMIT polynom. */
1192/* IERCOD: Error code, */
7fd59977 1193/* = 0, Ok */
0d969553
Y
1194/* = 1, required order of constraint is not managed here. */
1195/* COMMONS USED : */
7fd59977 1196/* ---------------- */
1197
0d969553 1198/* REFERENCES CALLED : */
7fd59977 1199/* ----------------------- */
1200
0d969553 1201/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 1202/* ----------------------------------- */
258ff83b 1203/* The part of HERMIT(*,2*i+j) table where j=1 or 2 and i=0 to IORDRE, */
0d969553
Y
1204/* contains the coefficients of the polynom of degree 2*IORDRE+1 */
1205/* such as ALL values in -1 and in +1 of this polynom and its */
1206/* derivatives till order of derivation IORDRE are NULL, */
1207/* EXCEPT for the derivative of order i: */
1208/* - valued 1 in -1 if j=1 */
1209/* - valued 1 in +1 if j=2. */
7fd59977 1210/* > */
1211/* **********************************************************************
1212*/
1213
0d969553 1214/* Name of the routine */
7fd59977 1215
1216
1217 /* Parameter adjustments */
1218 hermit_dim1 = (*iordre + 1) << 1;
1219 hermit_offset = hermit_dim1 + 1;
1220 hermit -= hermit_offset;
1221
1222 /* Function Body */
1223 ibb = AdvApp2Var_SysBase::mnfndeb_();
1224 if (ibb >= 3) {
1225 AdvApp2Var_SysBase::mgenmsg_("MMA1HER", 7L);
1226 }
1227 *iercod = 0;
1228
0d969553 1229/* --- Recover (IORDRE+2) coeff of 2*(IORDRE+1) Hermit polynoms --
7fd59977 1230*/
1231
1232 if (*iordre == 0) {
1233 hermit[hermit_dim1 + 1] = .5;
1234 hermit[hermit_dim1 + 2] = -.5;
1235
1236 hermit[(hermit_dim1 << 1) + 1] = .5;
1237 hermit[(hermit_dim1 << 1) + 2] = .5;
1238 } else if (*iordre == 1) {
1239 hermit[hermit_dim1 + 1] = .5;
1240 hermit[hermit_dim1 + 2] = -.75;
1241 hermit[hermit_dim1 + 3] = 0.;
1242 hermit[hermit_dim1 + 4] = .25;
1243
1244 hermit[(hermit_dim1 << 1) + 1] = .5;
1245 hermit[(hermit_dim1 << 1) + 2] = .75;
1246 hermit[(hermit_dim1 << 1) + 3] = 0.;
1247 hermit[(hermit_dim1 << 1) + 4] = -.25;
1248
1249 hermit[hermit_dim1 * 3 + 1] = .25;
1250 hermit[hermit_dim1 * 3 + 2] = -.25;
1251 hermit[hermit_dim1 * 3 + 3] = -.25;
1252 hermit[hermit_dim1 * 3 + 4] = .25;
1253
1254 hermit[(hermit_dim1 << 2) + 1] = -.25;
1255 hermit[(hermit_dim1 << 2) + 2] = -.25;
1256 hermit[(hermit_dim1 << 2) + 3] = .25;
1257 hermit[(hermit_dim1 << 2) + 4] = .25;
1258 } else if (*iordre == 2) {
1259 hermit[hermit_dim1 + 1] = .5;
1260 hermit[hermit_dim1 + 2] = -.9375;
1261 hermit[hermit_dim1 + 3] = 0.;
1262 hermit[hermit_dim1 + 4] = .625;
1263 hermit[hermit_dim1 + 5] = 0.;
1264 hermit[hermit_dim1 + 6] = -.1875;
1265
1266 hermit[(hermit_dim1 << 1) + 1] = .5;
1267 hermit[(hermit_dim1 << 1) + 2] = .9375;
1268 hermit[(hermit_dim1 << 1) + 3] = 0.;
1269 hermit[(hermit_dim1 << 1) + 4] = -.625;
1270 hermit[(hermit_dim1 << 1) + 5] = 0.;
1271 hermit[(hermit_dim1 << 1) + 6] = .1875;
1272
1273 hermit[hermit_dim1 * 3 + 1] = .3125;
1274 hermit[hermit_dim1 * 3 + 2] = -.4375;
1275 hermit[hermit_dim1 * 3 + 3] = -.375;
1276 hermit[hermit_dim1 * 3 + 4] = .625;
1277 hermit[hermit_dim1 * 3 + 5] = .0625;
1278 hermit[hermit_dim1 * 3 + 6] = -.1875;
1279
1280 hermit[(hermit_dim1 << 2) + 1] = -.3125;
1281 hermit[(hermit_dim1 << 2) + 2] = -.4375;
1282 hermit[(hermit_dim1 << 2) + 3] = .375;
1283 hermit[(hermit_dim1 << 2) + 4] = .625;
1284 hermit[(hermit_dim1 << 2) + 5] = -.0625;
1285 hermit[(hermit_dim1 << 2) + 6] = -.1875;
1286
1287 hermit[hermit_dim1 * 5 + 1] = .0625;
1288 hermit[hermit_dim1 * 5 + 2] = -.0625;
1289 hermit[hermit_dim1 * 5 + 3] = -.125;
1290 hermit[hermit_dim1 * 5 + 4] = .125;
1291 hermit[hermit_dim1 * 5 + 5] = .0625;
1292 hermit[hermit_dim1 * 5 + 6] = -.0625;
1293
1294 hermit[hermit_dim1 * 6 + 1] = .0625;
1295 hermit[hermit_dim1 * 6 + 2] = .0625;
1296 hermit[hermit_dim1 * 6 + 3] = -.125;
1297 hermit[hermit_dim1 * 6 + 4] = -.125;
1298 hermit[hermit_dim1 * 6 + 5] = .0625;
1299 hermit[hermit_dim1 * 6 + 6] = .0625;
1300 } else {
1301 *iercod = 1;
1302 }
1303
1304/* ------------------------------ The End -------------------------------
1305*/
1306
1307 AdvApp2Var_SysBase::maermsg_("MMA1HER", iercod, 7L);
1308 if (ibb >= 3) {
1309 AdvApp2Var_SysBase::mgsomsg_("MMA1HER", 7L);
1310 }
1311 return 0;
1312} /* mma1her_ */
1313//=======================================================================
1314//function : mma1jak_
1315//purpose :
1316//=======================================================================
1317int mma1jak_(integer *ndimen,
1318 integer *nbroot,
1319 integer *iordre,
1320 integer *ndgjac,
1321 doublereal *somtab,
1322 doublereal *diftab,
1323 doublereal *cgauss,
1324 doublereal *crvjac,
1325 integer *iercod)
1326{
1327 /* System generated locals */
1328 integer somtab_dim1, somtab_offset, diftab_dim1, diftab_offset,
96a95605 1329 crvjac_dim1, crvjac_offset;
41194117 1330
7fd59977 1331 /* Local variables */
1ef32e96 1332 integer ibb;
7fd59977 1333
1334/* **********************************************************************
1335*/
1336
0d969553 1337/* FUNCTION : */
7fd59977 1338/* ---------- */
0d969553
Y
1339/* Calculate the curve of approximation of a non-polynomial function */
1340/* in the base of Jacobi. */
7fd59977 1341
0d969553 1342/* KEYWORDS : */
7fd59977 1343/* ----------- */
0d969553 1344/* FUNCTION,DISCRETISATION,APPROXIMATION,CONSTRAINT,CURVE,JACOBI */
7fd59977 1345
0d969553 1346/* INPUT ARGUMENTS : */
7fd59977 1347/* ------------------ */
0d969553
Y
1348/* NDIMEN: Total dimension of the space (sum of dimensions */
1349/* of sub-spaces) */
258ff83b 1350/* NBROOT: Nb of points of discretization of the iso, extremities not */
0d969553
Y
1351/* included. */
1352/* IORDRE: Order of constraint at the extremities of the boundary */
1353/* -1 = no constraints, */
1354/* 0 = constraints of passage of limits (i.e. C0), */
1355/* 1 = C0 + constraints of 1st derivatives (i.e. C1), */
1356/* 2 = C1 + constraints of 2nd derivatives (i.e. C2). */
258ff83b 1357/* NDGJAC: Degree of development in series to be used for calculation in the */
0d969553
Y
1358/* base of Jacobi. */
1359
1360/* OUTPUT ARGUMENTS : */
7fd59977 1361/* ------------------- */
0d969553
Y
1362/* CRVJAC : Curve of approximation of FONCNP with (eventually) */
1363/* taking into account of constraints at the extremities. */
1364/* This curve is of degree NDGJAC. */
1365/* IERCOD : Error code : */
1366/* 0 = All is ok. */
1367/* 33 = Pb to return data of du block data */
1368/* of coeff. of integration by GAUSS method. */
1369/* by program MMAPPTT. */
1370
1371/* COMMONS USED : */
7fd59977 1372/* ---------------- */
1373
0d969553 1374/* REFERENCES CALLED : */
7fd59977 1375/* ----------------------- */
7fd59977 1376/* > */
1377/* **********************************************************************
1378*/
1379
0d969553 1380/* Name of the routine */
7fd59977 1381
1382 /* Parameter adjustments */
1383 diftab_dim1 = *nbroot / 2 + 1;
1384 diftab_offset = diftab_dim1;
1385 diftab -= diftab_offset;
1386 somtab_dim1 = *nbroot / 2 + 1;
1387 somtab_offset = somtab_dim1;
1388 somtab -= somtab_offset;
1389 crvjac_dim1 = *ndgjac + 1;
1390 crvjac_offset = crvjac_dim1;
1391 crvjac -= crvjac_offset;
7fd59977 1392
1393 /* Function Body */
1394 ibb = AdvApp2Var_SysBase::mnfndeb_();
1395 if (ibb >= 2) {
1396 AdvApp2Var_SysBase::mgenmsg_("MMA1JAK", 7L);
1397 }
1398 *iercod = 0;
1399
0d969553 1400/* ----------------- Recover coeffs of integration by Gauss -----------
7fd59977 1401*/
1402
1403 AdvApp2Var_ApproxF2var::mmapptt_(ndgjac, nbroot, iordre, cgauss, iercod);
1404 if (*iercod > 0) {
1405 *iercod = 33;
1406 goto L9999;
1407 }
1408
0d969553 1409/* --------------- Calculate the curve in the base of Jacobi -----------
7fd59977 1410*/
1411
1412 mmmapcoe_(ndimen, ndgjac, iordre, nbroot, &somtab[somtab_offset], &diftab[
1413 diftab_offset], cgauss, &crvjac[crvjac_offset]);
1414
1415/* ------------------------------ The End -------------------------------
1416*/
1417
1418L9999:
1419 if (*iercod != 0) {
1420 AdvApp2Var_SysBase::maermsg_("MMA1JAK", iercod, 7L);
1421 }
1422 if (ibb >= 2) {
1423 AdvApp2Var_SysBase::mgsomsg_("MMA1JAK", 7L);
1424 }
1425 return 0;
1426} /* mma1jak_ */
1427
1428//=======================================================================
1429//function : mma1noc_
1430//purpose :
1431//=======================================================================
1432int mma1noc_(doublereal *dfuvin,
1433 integer *ndimen,
1434 integer *iordre,
1435 doublereal *cntrin,
1436 doublereal *duvout,
1437 integer *isofav,
1438 integer *ideriv,
1439 doublereal *cntout)
1440{
1441 /* System generated locals */
1442 integer i__1;
1443 doublereal d__1;
41194117 1444
7fd59977 1445 /* Local variables */
1ef32e96
RL
1446 doublereal rider, riord;
1447 integer nd, ibb;
1448 doublereal bid;
7fd59977 1449/* **********************************************************************
1450*/
1451
0d969553 1452/* FUNCTION : */
7fd59977 1453/* ---------- */
0d969553
Y
1454/* Normalization of constraints of derivatives, defined on DFUVIN */
1455/* on block DUVOUT. */
7fd59977 1456
0d969553 1457/* KEYWORDS : */
7fd59977 1458/* ----------- */
0d969553 1459/* ALL, AB_SPECIFI::VECTEUR&,DERIVEE&,NORMALISATION,&VECTEUR */
7fd59977 1460
0d969553 1461/* INPUT ARGUMENTS : */
7fd59977 1462/* ------------------ */
0d969553 1463/* DFUVIN: Limits of the block of definition by U and by V where
7fd59977 1464*/
0d969553
Y
1465/* constraints CNTRIN are defined. */
1466/* NDIMEN: Dimension of the space. */
1467/* IORDRE: Order of constraint imposed at the extremities of the iso. */
1468/* (if Iso-U, it is necessary to calculate derivatives by V and vice */
7fd59977 1469/* versa). */
0d969553
Y
1470/* = 0, the extremities of the iso are calculated */
1471/* = 1, additionally the 1st derivative in the direction */
1472/* of the iso is calculated */
1473/* = 2, additionally the 2nd derivative in the direction */
1474/* of the iso is calculated */
1475/* CNTRIN: Contains, if IORDRE>=0, IORDRE+1 derivatives */
1476/* of order IORDRE of F(Uc,v) or of F(u,Vc), following the */
1477/* value of ISOFAV, RENORMALIZED by u and v in (-1,1). */
1478/* DUVOUT: Limits of the block of definition by U and by V where the */
1479/* constraints CNTOUT will be defined. */
1480/* ISOFAV: Isoparameter fixed for the discretization; */
1481/* = 1, discretization with fixed U=Uc and variable V. */
1482/* = 2, discretization with fixed V=Vc and variable U. */
7fd59977 1483/* IDERIV: Ordre de derivee transverse a l'iso fixee (Si Iso-U=Uc */
0d969553
Y
1484/* is fixed, the derivative of order IDERIV is discretized by U */
1485/* of F(Uc,v). The same if iso-V is fixed). */
1486/* Varies from (positioning) to 2 (2nd derivative). */
7fd59977 1487
0d969553 1488/* OUTPUT ARGUMENTS : */
7fd59977 1489/* ------------------- */
0d969553
Y
1490/* CNTOUT: Contains, if IORDRE>=0, IORDRE+1 derivatives */
1491/* of order IORDRE of F(Uc,v) or of F(u,Vc), depending on the */
1492/* value of ISOFAV, RENORMALIZED for u and v in DUVOUT. */
7fd59977 1493
0d969553 1494/* COMMONS USED : */
7fd59977 1495/* ---------------- */
1496
0d969553
Y
1497/* REFERENCES CALLED : */
1498/* --------------------- */
7fd59977 1499
0d969553
Y
1500/* DESCRIPTION/NOTES/LIMITATIONS : */
1501/* ------------------------------- */
1502/* CNTRIN can be an output/input argument, */
1503/* so the call: */
7fd59977 1504
1505/* CALL MMA1NOC(DFUVIN,NDIMEN,IORDRE,CNTRIN,DUVOUT */
1506/* 1 ,ISOFAV,IDERIV,CNTRIN) */
1507
0d969553 1508/* is correct. */
7fd59977 1509/* > */
1510/* **********************************************************************
1511*/
1512
0d969553 1513/* Name of the routine */
7fd59977 1514
1515
1516 /* Parameter adjustments */
1517 dfuvin -= 3;
1518 --cntout;
1519 --cntrin;
1520 duvout -= 3;
1521
1522 /* Function Body */
1523 ibb = AdvApp2Var_SysBase::mnfndeb_();
1524 if (ibb >= 3) {
1525 AdvApp2Var_SysBase::mgenmsg_("MMA1NOC", 7L);
1526 }
1527
0d969553 1528/* --------------- Determination of coefficients of normalization -------
7fd59977 1529 */
1530
1531 if (*isofav == 1) {
1532 d__1 = (dfuvin[4] - dfuvin[3]) / (duvout[4] - duvout[3]);
1533 rider = AdvApp2Var_MathBase::pow__di(&d__1, ideriv);
1534 d__1 = (dfuvin[6] - dfuvin[5]) / (duvout[6] - duvout[5]);
1535 riord = AdvApp2Var_MathBase::pow__di(&d__1, iordre);
1536
1537 } else {
1538 d__1 = (dfuvin[6] - dfuvin[5]) / (duvout[6] - duvout[5]);
1539 rider = AdvApp2Var_MathBase::pow__di(&d__1, ideriv);
1540 d__1 = (dfuvin[4] - dfuvin[3]) / (duvout[4] - duvout[3]);
1541 riord = AdvApp2Var_MathBase::pow__di(&d__1, iordre);
1542 }
1543
0d969553 1544/* ------------- Renormalization of the vector of constraint ---------------
7fd59977 1545*/
1546
1547 bid = rider * riord;
1548 i__1 = *ndimen;
1549 for (nd = 1; nd <= i__1; ++nd) {
1550 cntout[nd] = bid * cntrin[nd];
1551/* L100: */
1552 }
1553
1554/* ------------------------------ The end -------------------------------
1555*/
1556
1557 if (ibb >= 3) {
1558 AdvApp2Var_SysBase::mgsomsg_("MMA1NOC", 7L);
1559 }
1560 return 0;
1561} /* mma1noc_ */
1562
1563//=======================================================================
1564//function : mma1nop_
1565//purpose :
1566//=======================================================================
1567int mma1nop_(integer *nbroot,
1568 doublereal *rootlg,
1569 doublereal *uvfonc,
1570 integer *isofav,
1571 doublereal *ttable,
1572 integer *iercod)
1573
1574{
1575 /* System generated locals */
1576 integer i__1;
41194117 1577
7fd59977 1578 /* Local variables */
1ef32e96
RL
1579 doublereal alinu, blinu, alinv, blinv;
1580 integer ii, ibb;
7fd59977 1581
1582/* ***********************************************************************
1583 */
1584
0d969553 1585/* FUNCTION : */
7fd59977 1586/* ---------- */
0d969553
Y
1587/* Normalization of parameters of an iso, starting from */
1588/* parametric block and parameters on (-1,1). */
7fd59977 1589
0d969553 1590/* KEYWORDS : */
7fd59977 1591/* ----------- */
1592/* TOUS,AB_SPECIFI :: ISO&,POINT&,NORMALISATION,&POINT,&ISO */
1593
0d969553 1594/* INPUT ARGUMENTS : */
7fd59977 1595/* -------------------- */
0d969553
Y
1596/* NBROOT: Nb of points of discretisation INSIDE the iso */
1597/* defined on (-1,1). */
1598/* ROOTLG: Table of discretization parameters on )-1,1( */
1599/* of the iso. */
1600/* UVFONC: Block of definition of the iso */
1601/* ISOFAV: = 1, this is iso-u; =2, this is iso-v. */
1602
1603/* OUTPUT ARGUMENTS : */
7fd59977 1604/* --------------------- */
0d969553 1605/* TTABLE: Table of parameters renormalized on UVFONC of the iso.
7fd59977 1606*/
1607/* IERCOD: = 0, OK */
0d969553 1608/* = 1, ISOFAV is out of allowed values. */
7fd59977 1609
7fd59977 1610/* > */
1611/* **********************************************************************
1612*/
0d969553 1613/* Name of the routine */
7fd59977 1614
1615
1616 /* Parameter adjustments */
1617 --rootlg;
1618 uvfonc -= 3;
1619
1620 /* Function Body */
1621 ibb = AdvApp2Var_SysBase::mnfndeb_();
1622 if (ibb >= 3) {
1623 AdvApp2Var_SysBase::mgenmsg_("MMA1NOP", 7L);
1624 }
1625
1626 alinu = (uvfonc[4] - uvfonc[3]) / 2.;
1627 blinu = (uvfonc[4] + uvfonc[3]) / 2.;
1628 alinv = (uvfonc[6] - uvfonc[5]) / 2.;
1629 blinv = (uvfonc[6] + uvfonc[5]) / 2.;
1630
1631 if (*isofav == 1) {
1632 ttable[0] = uvfonc[5];
1633 i__1 = *nbroot;
1634 for (ii = 1; ii <= i__1; ++ii) {
1635 ttable[ii] = alinv * rootlg[ii] + blinv;
1636/* L100: */
1637 }
1638 ttable[*nbroot + 1] = uvfonc[6];
1639 } else if (*isofav == 2) {
1640 ttable[0] = uvfonc[3];
1641 i__1 = *nbroot;
1642 for (ii = 1; ii <= i__1; ++ii) {
1643 ttable[ii] = alinu * rootlg[ii] + blinu;
1644/* L200: */
1645 }
1646 ttable[*nbroot + 1] = uvfonc[4];
1647 } else {
1648 goto L9100;
1649 }
1650
1651 goto L9999;
1652
1653/* ------------------------------ THE END -------------------------------
1654*/
1655
1656L9100:
1657 *iercod = 1;
1658 goto L9999;
1659
1660L9999:
1661 if (*iercod != 0) {
1662 AdvApp2Var_SysBase::maermsg_("MMA1NOP", iercod, 7L);
1663 }
1664 if (ibb >= 3) {
1665 AdvApp2Var_SysBase::mgsomsg_("MMA1NOP", 7L);
1666 }
1667
1668 return 0 ;
1669
1670} /* mma1nop_ */
1671
1672//=======================================================================
1673//function : mma2ac1_
1674//purpose :
1675//=======================================================================
1676int AdvApp2Var_ApproxF2var::mma2ac1_(integer const *ndimen,
1677 integer const *mxujac,
1678 integer const *mxvjac,
1679 integer const *iordru,
1680 integer const *iordrv,
1681 doublereal const *contr1,
1682 doublereal const * contr2,
1683 doublereal const *contr3,
1684 doublereal const *contr4,
1685 doublereal const *uhermt,
1686 doublereal const *vhermt,
1687 doublereal *patjac)
1688
1689{
1690 /* System generated locals */
1691 integer contr1_dim1, contr1_dim2, contr1_offset, contr2_dim1, contr2_dim2,
1692 contr2_offset, contr3_dim1, contr3_dim2, contr3_offset,
1693 contr4_dim1, contr4_dim2, contr4_offset, uhermt_dim1,
1694 uhermt_offset, vhermt_dim1, vhermt_offset, patjac_dim1,
1695 patjac_dim2, patjac_offset, i__1, i__2, i__3, i__4, i__5;
41194117 1696
7fd59977 1697 /* Local variables */
1ef32e96
RL
1698 logical ldbg;
1699 integer ndgu, ndgv;
1700 doublereal bidu1, bidu2, bidv1, bidv2;
1701 integer ioru1, iorv1, ii, nd, jj, ku, kv;
1702 doublereal cnt1, cnt2, cnt3, cnt4;
7fd59977 1703
1704/* **********************************************************************
1705*/
1706
0d969553 1707/* FUNCTION : */
7fd59977 1708/* ---------- */
0d969553 1709/* Add polynoms of edge constraints. */
7fd59977 1710
0d969553 1711/* KEYWORDS : */
7fd59977 1712/* ----------- */
1713/* TOUS,AB_SPECIFI::POINT&,CONTRAINTE&,ADDITION,&POLYNOME */
1714
0d969553 1715/* INPUT ARGUMENTS : */
7fd59977 1716/* ------------------ */
0d969553
Y
1717/* NDIMEN: Dimension of the space. */
1718/* MXUJAC: Max degree of the polynom of approximation by U. The */
1719/* representation in the orthogonal base starts from degree */
1720/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
1721/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
1722/* MXVJAC: Max degree of the polynom of approximation by V. The */
1723/* representation in the orthogonal base starts from degree */
1724/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
1725/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
1726/* IORDRU: Order of the base of Jacobi (-1,0,1 or 2) by U. Corresponds */
1727/* to the step of constraints: C0, C1 or C2. */
1728/* IORDRV: Order of the base of Jacobi (-1,0,1 or 2) by V. Corresponds */
1729/* to the step of constraints: C0, C1 or C2. */
1730/* CONTR1: Contains, if IORDRU and IORDRV>=0, the values at the */
1731/* extremities of F(U0,V0) and its derivatives. */
1732/* CONTR2: Contains, if IORDRU and IORDRV>=0, the values at the */
1733/* extremities of F(U1,V0) and its derivatives. */
1734/* CONTR3: Contains, if IORDRU and IORDRV>=0, the values at the */
1735/* extremities of F(U0,V1) and its derivatives. */
1736/* CONTR4: Contains, if IORDRU and IORDRV>=0, the values at the */
1737/* extremities of F(U1,V1) and its derivatives. */
1738/* UHERMT: Coeff. of Hermit polynoms of order IORDRU. */
1739/* VHERMT: Coeff. of Hermit polynoms of order IORDRV. */
1740/* PATJAC: Table of coefficients of the polynom P(u,v) of approximation */
1741/* of F(u,v) WITHOUT taking into account the constraints. */
1742
1743/* OUTPUT ARGUMENTS : */
7fd59977 1744/* ------------------- */
0d969553
Y
1745/* PATJAC: Table of coefficients of the polynom P(u,v) by approximation */
1746/* of F(u,v) WITH taking into account of constraints. */
7fd59977 1747/* > */
1748/* **********************************************************************
1749*/
0d969553 1750/* Name of the routine */
7fd59977 1751
0d969553 1752/* --------------------------- Initialization --------------------------
7fd59977 1753*/
1754
1755 /* Parameter adjustments */
1756 patjac_dim1 = *mxujac + 1;
1757 patjac_dim2 = *mxvjac + 1;
1758 patjac_offset = patjac_dim1 * patjac_dim2;
1759 patjac -= patjac_offset;
1760 uhermt_dim1 = (*iordru << 1) + 2;
1761 uhermt_offset = uhermt_dim1;
1762 uhermt -= uhermt_offset;
1763 vhermt_dim1 = (*iordrv << 1) + 2;
1764 vhermt_offset = vhermt_dim1;
1765 vhermt -= vhermt_offset;
1766 contr4_dim1 = *ndimen;
1767 contr4_dim2 = *iordru + 2;
1768 contr4_offset = contr4_dim1 * (contr4_dim2 + 1) + 1;
1769 contr4 -= contr4_offset;
1770 contr3_dim1 = *ndimen;
1771 contr3_dim2 = *iordru + 2;
1772 contr3_offset = contr3_dim1 * (contr3_dim2 + 1) + 1;
1773 contr3 -= contr3_offset;
1774 contr2_dim1 = *ndimen;
1775 contr2_dim2 = *iordru + 2;
1776 contr2_offset = contr2_dim1 * (contr2_dim2 + 1) + 1;
1777 contr2 -= contr2_offset;
1778 contr1_dim1 = *ndimen;
1779 contr1_dim2 = *iordru + 2;
1780 contr1_offset = contr1_dim1 * (contr1_dim2 + 1) + 1;
1781 contr1 -= contr1_offset;
1782
1783 /* Function Body */
1784 ldbg = AdvApp2Var_SysBase::mnfndeb_() >= 3;
1785 if (ldbg) {
1786 AdvApp2Var_SysBase::mgenmsg_("MMA2AC1", 7L);
1787 }
1788
0d969553 1789/* ------------ SUBTRACTION OF ANGULAR CONSTRAINTS -------------------
7fd59977 1790*/
1791
1792 ioru1 = *iordru + 1;
1793 iorv1 = *iordrv + 1;
1794 ndgu = (*iordru << 1) + 1;
1795 ndgv = (*iordrv << 1) + 1;
1796
1797 i__1 = iorv1;
1798 for (jj = 1; jj <= i__1; ++jj) {
1799 i__2 = ioru1;
1800 for (ii = 1; ii <= i__2; ++ii) {
1801 i__3 = *ndimen;
1802 for (nd = 1; nd <= i__3; ++nd) {
1803 cnt1 = contr1[nd + (ii + jj * contr1_dim2) * contr1_dim1];
1804 cnt2 = contr2[nd + (ii + jj * contr2_dim2) * contr2_dim1];
1805 cnt3 = contr3[nd + (ii + jj * contr3_dim2) * contr3_dim1];
1806 cnt4 = contr4[nd + (ii + jj * contr4_dim2) * contr4_dim1];
1807 i__4 = ndgv;
1808 for (kv = 0; kv <= i__4; ++kv) {
1809 bidv1 = vhermt[kv + ((jj << 1) - 1) * vhermt_dim1];
1810 bidv2 = vhermt[kv + (jj << 1) * vhermt_dim1];
1811 i__5 = ndgu;
1812 for (ku = 0; ku <= i__5; ++ku) {
1813 bidu1 = uhermt[ku + ((ii << 1) - 1) * uhermt_dim1];
1814 bidu2 = uhermt[ku + (ii << 1) * uhermt_dim1];
1815 patjac[ku + (kv + nd * patjac_dim2) * patjac_dim1] =
1816 patjac[ku + (kv + nd * patjac_dim2) *
1817 patjac_dim1] - bidu1 * bidv1 * cnt1 - bidu2 *
1818 bidv1 * cnt2 - bidu1 * bidv2 * cnt3 - bidu2 *
1819 bidv2 * cnt4;
1820/* L500: */
1821 }
1822/* L400: */
1823 }
1824/* L300: */
1825 }
1826/* L200: */
1827 }
1828/* L100: */
1829 }
1830
1831/* ------------------------------ The end -------------------------------
1832*/
1833
1834 if (ldbg) {
1835 AdvApp2Var_SysBase::mgsomsg_("MMA2AC1", 7L);
1836 }
1837 return 0;
1838} /* mma2ac1_ */
1839
1840//=======================================================================
1841//function : mma2ac2_
1842//purpose :
1843//=======================================================================
1844int AdvApp2Var_ApproxF2var::mma2ac2_(const integer *ndimen,
1845 const integer *mxujac,
1846 const integer *mxvjac,
1847 const integer *iordrv,
1848 const integer *nclimu,
1849 const integer *ncfiv1,
1850 const doublereal *crbiv1,
1851 const integer *ncfiv2,
1852 const doublereal *crbiv2,
1853 const doublereal *vhermt,
1854 doublereal *patjac)
1855
1856{
1857 /* System generated locals */
1858 integer crbiv1_dim1, crbiv1_dim2, crbiv1_offset, crbiv2_dim1, crbiv2_dim2,
1859 crbiv2_offset, patjac_dim1, patjac_dim2, patjac_offset,
1860 vhermt_dim1, vhermt_offset, i__1, i__2, i__3, i__4;
41194117 1861
7fd59977 1862 /* Local variables */
1ef32e96
RL
1863 logical ldbg;
1864 integer ndgv1, ndgv2, ii, jj, nd, kk;
1865 doublereal bid1, bid2;
7fd59977 1866
1867/* **********************************************************************
1868*/
1869
0d969553 1870/* FUNCTION : */
7fd59977 1871/* ---------- */
0d969553 1872/* Add polynoms of constraints */
7fd59977 1873
0d969553 1874/* KEYWORDS : */
7fd59977 1875/* ----------- */
0d969553 1876/* FUNCTION,APPROXIMATION,COEFFICIENT,POLYNOM */
7fd59977 1877
0d969553 1878/* INPUT ARGUMENTS : */
7fd59977 1879/* ------------------ */
0d969553
Y
1880/* NDIMEN: Dimension of the space. */
1881/* MXUJAC: Max degree of the polynom of approximation by U. The */
1882/* representation in the orthogonal base starts from degree */
1883/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
1884/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
1885/* MXVJAC: Max degree of the polynom of approximation by V. The */
1886/* representation in the orthogonal base starts from degree */
1887/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
1888/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
1889/* IORDRV: Order of the base of Jacobi (-1,0,1 or 2) by V. Corresponds */
1890/* to the step of constraints: C0, C1 or C2. */
1891/* NCLIMU LIMIT nb of coeff by u of the solution P(u,v)
1892* NCFIV1: Nb of Coeff. of curves stored in CRBIV1. */
1893/* CRBIV1: Table of coeffs of the approximation of iso-V0 and its */
1894/* derivatives till order IORDRV. */
1895/* NCFIV2: Nb of Coeff. of curves stored in CRBIV2. */
1896/* CRBIV2: Table of coeffs of approximation of iso-V1 and its */
1897/* derivatives till order IORDRV. */
1898/* VHERMT: Coeff. of Hermit polynoms of order IORDRV. */
1899/* PATJAC: Table of coefficients of the polynom P(u,v) of approximation */
1900/* of F(u,v) WITHOUT taking into account the constraints. */
1901
1902/* OUTPUT ARGUMENTS : */
7fd59977 1903/* ------------------- */
0d969553
Y
1904/* PATJAC: Table of coefficients of the polynom P(u,v) by approximation */
1905/* of F(u,v) WITH taking into account of constraints. */
258ff83b 1906/* > */
7fd59977 1907
7fd59977 1908
7fd59977 1909/* > */
1910/* **********************************************************************
1911*/
0d969553 1912/* Name of the routine */
7fd59977 1913
1914/* --------------------------- Initialisations --------------------------
1915*/
1916
1917 /* Parameter adjustments */
1918 patjac_dim1 = *mxujac + 1;
1919 patjac_dim2 = *mxvjac + 1;
1920 patjac_offset = patjac_dim1 * patjac_dim2;
1921 patjac -= patjac_offset;
1922 vhermt_dim1 = (*iordrv << 1) + 2;
1923 vhermt_offset = vhermt_dim1;
1924 vhermt -= vhermt_offset;
1925 --ncfiv2;
1926 --ncfiv1;
1927 crbiv2_dim1 = *nclimu;
1928 crbiv2_dim2 = *ndimen;
1929 crbiv2_offset = crbiv2_dim1 * (crbiv2_dim2 + 1);
1930 crbiv2 -= crbiv2_offset;
1931 crbiv1_dim1 = *nclimu;
1932 crbiv1_dim2 = *ndimen;
1933 crbiv1_offset = crbiv1_dim1 * (crbiv1_dim2 + 1);
1934 crbiv1 -= crbiv1_offset;
1935
1936 /* Function Body */
1937 ldbg = AdvApp2Var_SysBase::mnfndeb_() >= 3;
1938 if (ldbg) {
1939 AdvApp2Var_SysBase::mgenmsg_("MMA2AC2", 7L);
1940 }
1941
0d969553 1942/* ------------ ADDING of coeff by u of curves, by v of Hermit --------
7fd59977 1943*/
1944
1945 i__1 = *iordrv + 1;
1946 for (ii = 1; ii <= i__1; ++ii) {
1947 ndgv1 = ncfiv1[ii] - 1;
1948 ndgv2 = ncfiv2[ii] - 1;
1949 i__2 = *ndimen;
1950 for (nd = 1; nd <= i__2; ++nd) {
1951 i__3 = (*iordrv << 1) + 1;
1952 for (jj = 0; jj <= i__3; ++jj) {
1953 bid1 = vhermt[jj + ((ii << 1) - 1) * vhermt_dim1];
1954 i__4 = ndgv1;
1955 for (kk = 0; kk <= i__4; ++kk) {
1956 patjac[kk + (jj + nd * patjac_dim2) * patjac_dim1] +=
1957 bid1 * crbiv1[kk + (nd + ii * crbiv1_dim2) *
1958 crbiv1_dim1];
1959/* L400: */
1960 }
1961 bid2 = vhermt[jj + (ii << 1) * vhermt_dim1];
1962 i__4 = ndgv2;
1963 for (kk = 0; kk <= i__4; ++kk) {
1964 patjac[kk + (jj + nd * patjac_dim2) * patjac_dim1] +=
1965 bid2 * crbiv2[kk + (nd + ii * crbiv2_dim2) *
1966 crbiv2_dim1];
1967/* L500: */
1968 }
1969/* L300: */
1970 }
1971/* L200: */
1972 }
1973/* L100: */
1974 }
1975
1976/* ------------------------------ The end -------------------------------
1977*/
1978
1979 if (ldbg) {
1980 AdvApp2Var_SysBase::mgsomsg_("MMA2AC2", 7L);
1981 }
1982 return 0;
1983} /* mma2ac2_ */
1984
1985
1986//=======================================================================
1987//function : mma2ac3_
1988//purpose :
1989//=======================================================================
1990int AdvApp2Var_ApproxF2var::mma2ac3_(const integer *ndimen,
1991 const integer *mxujac,
1992 const integer *mxvjac,
1993 const integer *iordru,
1994 const integer *nclimv,
1995 const integer *ncfiu1,
1996 const doublereal * crbiu1,
1997 const integer *ncfiu2,
1998 const doublereal *crbiu2,
1999 const doublereal *uhermt,
2000 doublereal *patjac)
2001
2002{
2003 /* System generated locals */
2004 integer crbiu1_dim1, crbiu1_dim2, crbiu1_offset, crbiu2_dim1, crbiu2_dim2,
2005 crbiu2_offset, patjac_dim1, patjac_dim2, patjac_offset,
2006 uhermt_dim1, uhermt_offset, i__1, i__2, i__3, i__4;
41194117 2007
7fd59977 2008 /* Local variables */
1ef32e96
RL
2009 logical ldbg;
2010 integer ndgu1, ndgu2, ii, jj, nd, kk;
2011 doublereal bid1, bid2;
7fd59977 2012
2013/* **********************************************************************
2014*/
2015
0d969553 2016/* FUNCTION : */
7fd59977 2017/* ---------- */
2018/* Ajout des polynomes de contraintes */
2019
0d969553 2020/* KEYWORDS : */
7fd59977 2021/* ----------- */
2022/* FONCTION,APPROXIMATION,COEFFICIENT,POLYNOME */
2023
0d969553 2024/* INPUT ARGUMENTS : */
7fd59977 2025/* ------------------ */
0d969553
Y
2026/* NDIMEN: Dimension of the space. */
2027/* MXUJAC: Max degree of the polynom of approximation by U. The */
2028/* representation in the orthogonal base starts from degree */
2029/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
2030/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
2031/* MXVJAC: Max degree of the polynom of approximation by V. The */
2032/* representation in the orthogonal base starts from degree */
2033/* 0 to degree MXUJAC-2*(IORDRU+1). The polynomial base is the */
2034/* base of Jacobi of order -1 (Legendre), 0, 1 or 2 */
2035/* IORDRU: Order of the base of Jacobi (-1,0,1 or 2) by U. Corresponds */
2036/* to the step of constraints: C0, C1 or C2. */
2037/* NCLIMV LIMIT nb of coeff by v of the solution P(u,v)
2038* NCFIU1: Nb of Coeff. of curves stored in CRBIU1. */
2039/* CRBIU1: Table of coeffs of the approximation of iso-U0 and its */
2040/* derivatives till order IORDRU. */
2041/* NCFIU2: Nb of Coeff. of curves stored in CRBIU2. */
2042/* CRBIU2: Table of coeffs of approximation of iso-U1 and its */
2043/* derivatives till order IORDRU */
2044/* UHERMT: Coeff. of Hermit polynoms of order IORDRU. */
2045/* PATJAC: Table of coefficients of the polynom P(u,v) of approximation */
2046/* of F(u,v) WITHOUT taking into account the constraints. */
2047
2048/* OUTPUT ARGUMENTS : */
7fd59977 2049/* ------------------- */
0d969553
Y
2050/* PATJAC: Table of coefficients of the polynom P(u,v) by approximation */
2051/* of F(u,v) WITH taking into account of constraints. */
7fd59977 2052
7fd59977 2053
7fd59977 2054/* > */
2055/* **********************************************************************
2056*/
0d969553 2057/* The name of the routine */
7fd59977 2058
0d969553 2059/* --------------------------- Initializations --------------------------
7fd59977 2060*/
2061
2062 /* Parameter adjustments */
2063 patjac_dim1 = *mxujac + 1;
2064 patjac_dim2 = *mxvjac + 1;
2065 patjac_offset = patjac_dim1 * patjac_dim2;
2066 patjac -= patjac_offset;
2067 uhermt_dim1 = (*iordru << 1) + 2;
2068 uhermt_offset = uhermt_dim1;
2069 uhermt -= uhermt_offset;
2070 --ncfiu2;
2071 --ncfiu1;
2072 crbiu2_dim1 = *nclimv;
2073 crbiu2_dim2 = *ndimen;
2074 crbiu2_offset = crbiu2_dim1 * (crbiu2_dim2 + 1);
2075 crbiu2 -= crbiu2_offset;
2076 crbiu1_dim1 = *nclimv;
2077 crbiu1_dim2 = *ndimen;
2078 crbiu1_offset = crbiu1_dim1 * (crbiu1_dim2 + 1);
2079 crbiu1 -= crbiu1_offset;
2080
2081 /* Function Body */
2082 ldbg = AdvApp2Var_SysBase::mnfndeb_() >= 3;
2083 if (ldbg) {
2084 AdvApp2Var_SysBase::mgenmsg_("MMA2AC3", 7L);
2085 }
2086
0d969553 2087/* ------------ ADDING of coeff by u of curves, by v of Hermit --------
7fd59977 2088*/
2089
2090 i__1 = *iordru + 1;
2091 for (ii = 1; ii <= i__1; ++ii) {
2092 ndgu1 = ncfiu1[ii] - 1;
2093 ndgu2 = ncfiu2[ii] - 1;
2094 i__2 = *ndimen;
2095 for (nd = 1; nd <= i__2; ++nd) {
2096 i__3 = ndgu1;
2097 for (jj = 0; jj <= i__3; ++jj) {
2098 bid1 = crbiu1[jj + (nd + ii * crbiu1_dim2) * crbiu1_dim1];
2099 i__4 = (*iordru << 1) + 1;
2100 for (kk = 0; kk <= i__4; ++kk) {
2101 patjac[kk + (jj + nd * patjac_dim2) * patjac_dim1] +=
2102 bid1 * uhermt[kk + ((ii << 1) - 1) * uhermt_dim1];
2103/* L400: */
2104 }
2105/* L300: */
2106 }
2107 i__3 = ndgu2;
2108 for (jj = 0; jj <= i__3; ++jj) {
2109 bid2 = crbiu2[jj + (nd + ii * crbiu2_dim2) * crbiu2_dim1];
2110 i__4 = (*iordru << 1) + 1;
2111 for (kk = 0; kk <= i__4; ++kk) {
2112 patjac[kk + (jj + nd * patjac_dim2) * patjac_dim1] +=
2113 bid2 * uhermt[kk + (ii << 1) * uhermt_dim1];
2114/* L600: */
2115 }
2116/* L500: */
2117 }
2118
2119/* L200: */
2120 }
2121/* L100: */
2122 }
2123
2124/* ------------------------------ The end -------------------------------
2125*/
2126
2127 if (ldbg) {
2128 AdvApp2Var_SysBase::mgsomsg_("MMA2AC3", 7L);
2129 }
2130 return 0;
2131} /* mma2ac3_ */
2132
2133//=======================================================================
2134//function : mma2can_
2135//purpose :
2136//=======================================================================
2137int AdvApp2Var_ApproxF2var::mma2can_(const integer *ncfmxu,
2138 const integer *ncfmxv,
2139 const integer *ndimen,
2140 const integer *iordru,
2141 const integer *iordrv,
2142 const integer *ncoefu,
2143 const integer *ncoefv,
2144 const doublereal *patjac,
2145 doublereal *pataux,
2146 doublereal *patcan,
2147 integer *iercod)
2148
2149{
2150 /* System generated locals */
2151 integer patjac_dim1, patjac_dim2, patjac_offset, patcan_dim1, patcan_dim2,
2152 patcan_offset, i__1, i__2;
41194117 2153
7fd59977 2154 /* Local variables */
1ef32e96
RL
2155 logical ldbg;
2156 integer ilon1, ilon2, ii, nd;
7fd59977 2157
2158/* **********************************************************************
2159*/
2160
0d969553 2161/* FUNCTION : */
7fd59977 2162/* ---------- */
0d969553
Y
2163/* Change of Jacobi base to canonical (-1,1) and writing in a greater */
2164/* table. */
7fd59977 2165
0d969553 2166/* KEYWORDS : */
7fd59977 2167/* ----------- */
0d969553 2168/* ALL,AB_SPECIFI,CARREAU&,CONVERSION,JACOBI,CANNONIQUE,&CARREAU */
7fd59977 2169
0d969553 2170/* INPUT ARGUMENTS : */
7fd59977 2171/* -------------------- */
0d969553
Y
2172/* NCFMXU: Dimension by U of resulting table PATCAN */
2173/* NCFMXV: Dimension by V of resulting table PATCAN */
2174/* NDIMEN: Dimension of the workspace. */
2175/* IORDRU: Order of constraint by U */
2176/* IORDRV: Order of constraint by V. */
2177/* NCOEFU: Nb of coeff by U of square PATJAC */
2178/* NCOEFV: Nb of coeff by V of square PATJAC */
2179/* PATJAC: Square in the base of Jacobi of order IORDRU by U and */
2180/* IORDRV by V. */
2181
2182/* OUTPUT ARGUMENTS : */
7fd59977 2183/* --------------------- */
0d969553
Y
2184/* PATAUX: Auxiliary Table. */
2185/* PATCAN: Table of coefficients in the canonic base. */
2186/* IERCOD: Error code. */
2187/* = 0, everything goes well, and all things are equal. */
2188/* = 1, the program refuses to process with incorrect input arguments */
2189
2190
2191/* COMMONS USED : */
7fd59977 2192/* ------------------ */
2193
0d969553 2194/* REFERENCES CALLED : */
7fd59977 2195/* --------------------- */
2196
0d969553 2197/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 2198/* ----------------------------------- */
7fd59977 2199/* > */
2200/* **********************************************************************
2201*/
2202
2203
2204 /* Parameter adjustments */
2205 patcan_dim1 = *ncfmxu;
2206 patcan_dim2 = *ncfmxv;
2207 patcan_offset = patcan_dim1 * (patcan_dim2 + 1) + 1;
2208 patcan -= patcan_offset;
2209 --pataux;
2210 patjac_dim1 = *ncoefu;
2211 patjac_dim2 = *ncoefv;
2212 patjac_offset = patjac_dim1 * (patjac_dim2 + 1) + 1;
2213 patjac -= patjac_offset;
2214
2215 /* Function Body */
2216 ldbg = AdvApp2Var_SysBase::mnfndeb_() >= 2;
2217 if (ldbg) {
2218 AdvApp2Var_SysBase::mgenmsg_("MMA2CAN", 7L);
2219 }
2220 *iercod = 0;
2221
2222 if (*iordru < -1 || *iordru > 2) {
2223 goto L9100;
2224 }
2225 if (*iordrv < -1 || *iordrv > 2) {
2226 goto L9100;
2227 }
2228 if (*ncoefu > *ncfmxu || *ncoefv > *ncfmxv) {
2229 goto L9100;
2230 }
2231
0d969553 2232/* --> Pass to canonic base (-1,1) */
7fd59977 2233 mmjacpt_(ndimen, ncoefu, ncoefv, iordru, iordrv, &patjac[patjac_offset], &
2234 pataux[1], &patcan[patcan_offset]);
2235
0d969553 2236/* --> Write all in a greater table */
fadcea2c
RL
2237 AdvApp2Var_MathBase::mmfmca8_(ncoefu,
2238 ncoefv,
2239 ndimen,
2240 ncfmxu,
2241 ncfmxv,
2242 ndimen,
2243 &patcan[patcan_offset],
2244 &patcan[patcan_offset]);
7fd59977 2245
0d969553 2246/* --> Complete with zeros the resulting table. */
7fd59977 2247 ilon1 = *ncfmxu - *ncoefu;
2248 ilon2 = *ncfmxu * (*ncfmxv - *ncoefv);
2249 i__1 = *ndimen;
2250 for (nd = 1; nd <= i__1; ++nd) {
2251 if (ilon1 > 0) {
2252 i__2 = *ncoefv;
2253 for (ii = 1; ii <= i__2; ++ii) {
2254 AdvApp2Var_SysBase::mvriraz_(&ilon1,
fadcea2c 2255 &patcan[*ncoefu + 1 + (ii + nd * patcan_dim2) * patcan_dim1]);
7fd59977 2256/* L110: */
2257 }
2258 }
2259 if (ilon2 > 0) {
2260 AdvApp2Var_SysBase::mvriraz_(&ilon2,
fadcea2c 2261 &patcan[(*ncoefv + 1 + nd * patcan_dim2) * patcan_dim1 + 1]);
7fd59977 2262 }
2263/* L100: */
2264 }
2265
2266 goto L9999;
2267
0d969553 2268/* ----------------------
7fd59977 2269*/
2270
2271L9100:
2272 *iercod = 1;
2273 goto L9999;
2274
2275L9999:
2276 AdvApp2Var_SysBase::maermsg_("MMA2CAN", iercod, 7L);
2277 if (ldbg) {
2278 AdvApp2Var_SysBase::mgsomsg_("MMA2CAN", 7L);
2279 }
2280 return 0 ;
2281} /* mma2can_ */
2282
2283//=======================================================================
2284//function : mma2cd1_
2285//purpose :
2286//=======================================================================
2287int mma2cd1_(integer *ndimen,
2288 integer *nbpntu,
2289 doublereal *urootl,
2290 integer *nbpntv,
2291 doublereal *vrootl,
2292 integer *iordru,
2293 integer *iordrv,
2294 doublereal *contr1,
2295 doublereal *contr2,
2296 doublereal *contr3,
2297 doublereal *contr4,
2298 doublereal *fpntbu,
2299 doublereal *fpntbv,
2300 doublereal *uhermt,
2301 doublereal *vhermt,
2302 doublereal *sosotb,
2303 doublereal *soditb,
2304 doublereal *disotb,
2305 doublereal *diditb)
2306
2307{
1ef32e96 2308 integer c__1 = 1;
41194117 2309
7fd59977 2310/* System generated locals */
2311 integer contr1_dim1, contr1_dim2, contr1_offset, contr2_dim1, contr2_dim2,
2312 contr2_offset, contr3_dim1, contr3_dim2, contr3_offset,
2313 contr4_dim1, contr4_dim2, contr4_offset, uhermt_dim1,
2314 uhermt_offset, vhermt_dim1, vhermt_offset, fpntbu_dim1,
2315 fpntbu_offset, fpntbv_dim1, fpntbv_offset, sosotb_dim1,
2316 sosotb_dim2, sosotb_offset, diditb_dim1, diditb_dim2,
2317 diditb_offset, soditb_dim1, soditb_dim2, soditb_offset,
2318 disotb_dim1, disotb_dim2, disotb_offset, i__1, i__2, i__3, i__4,
2319 i__5;
2320
2321 /* Local variables */
1ef32e96 2322 integer ncfhu, ncfhv, nuroo, nvroo, nd, ii, jj, kk, ll, ibb, kkm,
7fd59977 2323 llm, kkp, llp;
1ef32e96
RL
2324 doublereal bid1, bid2, bid3, bid4;
2325 doublereal diu1, diu2, div1, div2, sou1, sou2, sov1, sov2;
7fd59977 2326
7fd59977 2327/* **********************************************************************
2328*/
2329
0d969553 2330/* FUNCTION : */
7fd59977 2331/* ---------- */
0d969553
Y
2332/* Discretisation on the parameters of polynoms of interpolation */
2333/* of constraints at the corners of order IORDRE. */
7fd59977 2334
0d969553 2335/* KEYWORDS : */
7fd59977 2336/* ----------- */
2337/* TOUS, AB_SPECIFI::CONTRAINTE&, DISCRETISATION, &POINT */
2338
0d969553 2339/* INPUT ARGUMENTS : */
7fd59977 2340/* ------------------ */
0d969553
Y
2341/* NDIMEN: Dimension of the space. */
2342/* NBPNTU: Nb of INTERNAL parameters of discretisation by U. */
2343/* This is also the nb of root of Legendre polynom where discretization is done. */
2344/* UROOTL: Table of parameters of discretisation ON (-1,1) by U.
2345*/
2346/* NBPNTV: Nb of INTERNAL parameters of discretisation by V. */
2347/* This is also the nb of root of Legendre polynom where discretization is done. */
258ff83b 2348/* VROOTL: Table of discretization parameters on (-1,1) by V. */
0d969553
Y
2349/* IORDRU: Order of constraint imposed at the extremities of iso-V */
2350/* = 0, calculate the extremities of iso-V */
2351/* = 1, calculate, additionally, the 1st derivative in the direction of iso-V */
2352/* = 2, calculate, additionally, the 2nd derivative in the direction of iso-V */
2353/* IORDRV: Order of constraint imposed at the extremities of iso-U */
2354/* = 0, calculate the extremities of iso-U */
2355/* = 1, calculate, additionally, the 1st derivative in the direction of iso-U */
2356/* = 2, calculate, additionally, the 2nd derivative in the direction of iso-U */
2357/* CONTR1: Contains, if IORDRU and IORDRV>=0, the values at the */
2358/* extremities of F(U0,V0) and its derivatives. */
2359/* CONTR2: Contains, if IORDRU and IORDRV>=0, the values at the */
2360/* extremities of F(U1,V0) and its derivatives. */
2361/* CONTR3: Contains, if IORDRU and IORDRV>=0, the values at the */
2362/* extremities of F(U0,V1) and its derivatives. */
2363/* CONTR4: Contains, if IORDRU and IORDRV>=0, the values at the */
2364/* extremities of F(U1,V1) and its derivatives. */
2365/* SOSOTB: Preinitialized table (input/output argument). */
2366/* DISOTB: Preinitialized table (input/output argument). */
2367/* SODITB: Preinitialized table (input/output argument). */
2368/* DIDITB: Preinitialized table (input/output argument) */
2369
2370/* OUTPUT ARGUMENTS : */
7fd59977 2371/* ------------------- */
0d969553
Y
2372/* FPNTBU: Auxiliary table. */
2373/* FPNTBV: Auxiliary table. */
2374/* UHERMT: Table of 2*(IORDRU+1) coeff. of 2*(IORDRU+1) polynoms of Hermite. */
2375/* VHERMT: Table of 2*(IORDRV+1) coeff. of 2*(IORDRV+1) polynoms of Hermite. */
2376/* SOSOTB: Table where the terms of constraints are added */
7fd59977 2377/* C(ui,vj) + C(ui,-vj) + C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
2378/* with ui and vj positive roots of the Legendre polynom */
2379/* of degree NBPNTU and NBPNTV respectively. */
2380/* DISOTB: Table where the terms of constraints are added */
7fd59977 2381/* C(ui,vj) + C(ui,-vj) - C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
2382/* with ui and vj positive roots of the polynom of Legendre */
2383/* of degree NBPNTU and NBPNTV respectively. */
2384/* SODITB: Table where the terms of constraints are added */
7fd59977 2385/* C(ui,vj) - C(ui,-vj) + C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
2386/* with ui and vj positive roots of the polynom of Legendre */
2387/* of degree NBPNTU and NBPNTV respectively. */
2388/* DIDITB: Table where the terms of constraints are added */
7fd59977 2389/* C(ui,vj) - C(ui,-vj) - C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
2390/* with ui and vj positive roots of the polynom of Legendre */
2391/* of degree NBPNTU and NBPNTV respectively. */
7fd59977 2392
0d969553 2393/* COMMONS USED : */
7fd59977 2394/* ---------------- */
2395
0d969553 2396/* REFERENCES CALLED : */
7fd59977 2397/* ----------------------- */
2398
0d969553 2399/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 2400/* ----------------------------------- */
2401
7fd59977 2402/* > */
2403/* **********************************************************************
2404*/
2405
0d969553 2406/* Name of the routine */
7fd59977 2407
2408
2409 /* Parameter adjustments */
2410 --urootl;
2411 diditb_dim1 = *nbpntu / 2 + 1;
2412 diditb_dim2 = *nbpntv / 2 + 1;
2413 diditb_offset = diditb_dim1 * diditb_dim2;
2414 diditb -= diditb_offset;
2415 disotb_dim1 = *nbpntu / 2;
2416 disotb_dim2 = *nbpntv / 2;
2417 disotb_offset = disotb_dim1 * (disotb_dim2 + 1) + 1;
2418 disotb -= disotb_offset;
2419 soditb_dim1 = *nbpntu / 2;
2420 soditb_dim2 = *nbpntv / 2;
2421 soditb_offset = soditb_dim1 * (soditb_dim2 + 1) + 1;
2422 soditb -= soditb_offset;
2423 sosotb_dim1 = *nbpntu / 2 + 1;
2424 sosotb_dim2 = *nbpntv / 2 + 1;
2425 sosotb_offset = sosotb_dim1 * sosotb_dim2;
2426 sosotb -= sosotb_offset;
2427 --vrootl;
2428 uhermt_dim1 = (*iordru << 1) + 2;
2429 uhermt_offset = uhermt_dim1;
2430 uhermt -= uhermt_offset;
2431 fpntbu_dim1 = *nbpntu;
2432 fpntbu_offset = fpntbu_dim1 + 1;
2433 fpntbu -= fpntbu_offset;
2434 vhermt_dim1 = (*iordrv << 1) + 2;
2435 vhermt_offset = vhermt_dim1;
2436 vhermt -= vhermt_offset;
2437 fpntbv_dim1 = *nbpntv;
2438 fpntbv_offset = fpntbv_dim1 + 1;
2439 fpntbv -= fpntbv_offset;
2440 contr4_dim1 = *ndimen;
2441 contr4_dim2 = *iordru + 2;
2442 contr4_offset = contr4_dim1 * (contr4_dim2 + 1) + 1;
2443 contr4 -= contr4_offset;
2444 contr3_dim1 = *ndimen;
2445 contr3_dim2 = *iordru + 2;
2446 contr3_offset = contr3_dim1 * (contr3_dim2 + 1) + 1;
2447 contr3 -= contr3_offset;
2448 contr2_dim1 = *ndimen;
2449 contr2_dim2 = *iordru + 2;
2450 contr2_offset = contr2_dim1 * (contr2_dim2 + 1) + 1;
2451 contr2 -= contr2_offset;
2452 contr1_dim1 = *ndimen;
2453 contr1_dim2 = *iordru + 2;
2454 contr1_offset = contr1_dim1 * (contr1_dim2 + 1) + 1;
2455 contr1 -= contr1_offset;
2456
2457 /* Function Body */
2458 ibb = AdvApp2Var_SysBase::mnfndeb_();
2459 if (ibb >= 3) {
2460 AdvApp2Var_SysBase::mgenmsg_("MMA2CD1", 7L);
2461 }
2462
0d969553 2463/* ------------------- Discretisation of Hermite polynoms -----------
7fd59977 2464*/
2465
2466 ncfhu = (*iordru + 1) << 1;
2467 i__1 = ncfhu;
2468 for (ii = 1; ii <= i__1; ++ii) {
2469 i__2 = *nbpntu;
2470 for (ll = 1; ll <= i__2; ++ll) {
2471 AdvApp2Var_MathBase::mmmpocur_(&ncfhu, &c__1, &ncfhu, &uhermt[ii * uhermt_dim1], &
2472 urootl[ll], &fpntbu[ll + ii * fpntbu_dim1]);
2473/* L20: */
2474 }
2475/* L10: */
2476 }
2477 ncfhv = (*iordrv + 1) << 1;
2478 i__1 = ncfhv;
2479 for (jj = 1; jj <= i__1; ++jj) {
2480 i__2 = *nbpntv;
2481 for (kk = 1; kk <= i__2; ++kk) {
2482 AdvApp2Var_MathBase::mmmpocur_(&ncfhv, &c__1, &ncfhv, &vhermt[jj * vhermt_dim1], &
2483 vrootl[kk], &fpntbv[kk + jj * fpntbv_dim1]);
2484/* L40: */
2485 }
2486/* L30: */
2487 }
2488
0d969553 2489/* ---- The discretizations of polynoms of constraints are subtracted ----
7fd59977 2490*/
2491
2492 nuroo = *nbpntu / 2;
2493 nvroo = *nbpntv / 2;
2494 i__1 = *ndimen;
2495 for (nd = 1; nd <= i__1; ++nd) {
2496
2497 i__2 = *iordrv + 1;
2498 for (jj = 1; jj <= i__2; ++jj) {
2499 i__3 = *iordru + 1;
2500 for (ii = 1; ii <= i__3; ++ii) {
2501 bid1 = contr1[nd + (ii + jj * contr1_dim2) * contr1_dim1];
2502 bid2 = contr2[nd + (ii + jj * contr2_dim2) * contr2_dim1];
2503 bid3 = contr3[nd + (ii + jj * contr3_dim2) * contr3_dim1];
2504 bid4 = contr4[nd + (ii + jj * contr4_dim2) * contr4_dim1];
2505
2506 i__4 = nvroo;
2507 for (kk = 1; kk <= i__4; ++kk) {
2508 kkp = (*nbpntv + 1) / 2 + kk;
2509 kkm = nvroo - kk + 1;
2510 sov1 = fpntbv[kkp + ((jj << 1) - 1) * fpntbv_dim1] +
2511 fpntbv[kkm + ((jj << 1) - 1) * fpntbv_dim1];
2512 div1 = fpntbv[kkp + ((jj << 1) - 1) * fpntbv_dim1] -
2513 fpntbv[kkm + ((jj << 1) - 1) * fpntbv_dim1];
2514 sov2 = fpntbv[kkp + (jj << 1) * fpntbv_dim1] + fpntbv[kkm
2515 + (jj << 1) * fpntbv_dim1];
2516 div2 = fpntbv[kkp + (jj << 1) * fpntbv_dim1] - fpntbv[kkm
2517 + (jj << 1) * fpntbv_dim1];
2518 i__5 = nuroo;
2519 for (ll = 1; ll <= i__5; ++ll) {
2520 llp = (*nbpntu + 1) / 2 + ll;
2521 llm = nuroo - ll + 1;
2522 sou1 = fpntbu[llp + ((ii << 1) - 1) * fpntbu_dim1] +
2523 fpntbu[llm + ((ii << 1) - 1) * fpntbu_dim1];
2524 diu1 = fpntbu[llp + ((ii << 1) - 1) * fpntbu_dim1] -
2525 fpntbu[llm + ((ii << 1) - 1) * fpntbu_dim1];
2526 sou2 = fpntbu[llp + (ii << 1) * fpntbu_dim1] + fpntbu[
2527 llm + (ii << 1) * fpntbu_dim1];
2528 diu2 = fpntbu[llp + (ii << 1) * fpntbu_dim1] - fpntbu[
2529 llm + (ii << 1) * fpntbu_dim1];
2530 sosotb[ll + (kk + nd * sosotb_dim2) * sosotb_dim1] =
2531 sosotb[ll + (kk + nd * sosotb_dim2) *
2532 sosotb_dim1] - bid1 * sou1 * sov1 - bid2 *
2533 sou2 * sov1 - bid3 * sou1 * sov2 - bid4 *
2534 sou2 * sov2;
2535 soditb[ll + (kk + nd * soditb_dim2) * soditb_dim1] =
2536 soditb[ll + (kk + nd * soditb_dim2) *
2537 soditb_dim1] - bid1 * sou1 * div1 - bid2 *
2538 sou2 * div1 - bid3 * sou1 * div2 - bid4 *
2539 sou2 * div2;
2540 disotb[ll + (kk + nd * disotb_dim2) * disotb_dim1] =
2541 disotb[ll + (kk + nd * disotb_dim2) *
2542 disotb_dim1] - bid1 * diu1 * sov1 - bid2 *
2543 diu2 * sov1 - bid3 * diu1 * sov2 - bid4 *
2544 diu2 * sov2;
2545 diditb[ll + (kk + nd * diditb_dim2) * diditb_dim1] =
2546 diditb[ll + (kk + nd * diditb_dim2) *
2547 diditb_dim1] - bid1 * diu1 * div1 - bid2 *
2548 diu2 * div1 - bid3 * diu1 * div2 - bid4 *
2549 diu2 * div2;
2550/* L450: */
2551 }
2552/* L400: */
2553 }
2554
0d969553
Y
2555/* ------------ Case when the discretization is done only on the roots
2556----------- */
2557/* ---------- of Legendre polynom of uneven degree, 0 is root
7fd59977 2558----------- */
7fd59977 2559
2560 if (*nbpntu % 2 == 1) {
2561 sou1 = fpntbu[nuroo + 1 + ((ii << 1) - 1) * fpntbu_dim1];
2562 sou2 = fpntbu[nuroo + 1 + (ii << 1) * fpntbu_dim1];
2563 i__4 = nvroo;
2564 for (kk = 1; kk <= i__4; ++kk) {
2565 kkp = (*nbpntv + 1) / 2 + kk;
2566 kkm = nvroo - kk + 1;
2567 sov1 = fpntbv[kkp + ((jj << 1) - 1) * fpntbv_dim1] +
2568 fpntbv[kkm + ((jj << 1) - 1) * fpntbv_dim1];
2569 div1 = fpntbv[kkp + ((jj << 1) - 1) * fpntbv_dim1] -
2570 fpntbv[kkm + ((jj << 1) - 1) * fpntbv_dim1];
2571 sov2 = fpntbv[kkp + (jj << 1) * fpntbv_dim1] + fpntbv[
2572 kkm + (jj << 1) * fpntbv_dim1];
2573 div2 = fpntbv[kkp + (jj << 1) * fpntbv_dim1] - fpntbv[
2574 kkm + (jj << 1) * fpntbv_dim1];
2575 sosotb[(kk + nd * sosotb_dim2) * sosotb_dim1] =
2576 sosotb[(kk + nd * sosotb_dim2) * sosotb_dim1]
2577 - bid1 * sou1 * sov1 - bid2 * sou2 * sov1 -
2578 bid3 * sou1 * sov2 - bid4 * sou2 * sov2;
2579 diditb[(kk + nd * diditb_dim2) * diditb_dim1] =
2580 diditb[(kk + nd * diditb_dim2) * diditb_dim1]
2581 - bid1 * sou1 * div1 - bid2 * sou2 * div1 -
2582 bid3 * sou1 * div2 - bid4 * sou2 * div2;
2583/* L500: */
2584 }
2585 }
2586
2587 if (*nbpntv % 2 == 1) {
2588 sov1 = fpntbv[nvroo + 1 + ((jj << 1) - 1) * fpntbv_dim1];
2589 sov2 = fpntbv[nvroo + 1 + (jj << 1) * fpntbv_dim1];
2590 i__4 = nuroo;
2591 for (ll = 1; ll <= i__4; ++ll) {
2592 llp = (*nbpntu + 1) / 2 + ll;
2593 llm = nuroo - ll + 1;
2594 sou1 = fpntbu[llp + ((ii << 1) - 1) * fpntbu_dim1] +
2595 fpntbu[llm + ((ii << 1) - 1) * fpntbu_dim1];
2596 diu1 = fpntbu[llp + ((ii << 1) - 1) * fpntbu_dim1] -
2597 fpntbu[llm + ((ii << 1) - 1) * fpntbu_dim1];
2598 sou2 = fpntbu[llp + (ii << 1) * fpntbu_dim1] + fpntbu[
2599 llm + (ii << 1) * fpntbu_dim1];
2600 diu2 = fpntbu[llp + (ii << 1) * fpntbu_dim1] - fpntbu[
2601 llm + (ii << 1) * fpntbu_dim1];
2602 sosotb[ll + nd * sosotb_dim2 * sosotb_dim1] = sosotb[
2603 ll + nd * sosotb_dim2 * sosotb_dim1] - bid1 *
2604 sou1 * sov1 - bid2 * sou2 * sov1 - bid3 *
2605 sou1 * sov2 - bid4 * sou2 * sov2;
2606 diditb[ll + nd * diditb_dim2 * diditb_dim1] = diditb[
2607 ll + nd * diditb_dim2 * diditb_dim1] - bid1 *
2608 diu1 * sov1 - bid2 * diu2 * sov1 - bid3 *
2609 diu1 * sov2 - bid4 * diu2 * sov2;
2610/* L600: */
2611 }
2612 }
2613
2614 if (*nbpntu % 2 == 1 && *nbpntv % 2 == 1) {
2615 sou1 = fpntbu[nuroo + 1 + ((ii << 1) - 1) * fpntbu_dim1];
2616 sou2 = fpntbu[nuroo + 1 + (ii << 1) * fpntbu_dim1];
2617 sov1 = fpntbv[nvroo + 1 + ((jj << 1) - 1) * fpntbv_dim1];
2618 sov2 = fpntbv[nvroo + 1 + (jj << 1) * fpntbv_dim1];
2619 sosotb[nd * sosotb_dim2 * sosotb_dim1] = sosotb[nd *
2620 sosotb_dim2 * sosotb_dim1] - bid1 * sou1 * sov1 -
2621 bid2 * sou2 * sov1 - bid3 * sou1 * sov2 - bid4 *
2622 sou2 * sov2;
2623 diditb[nd * diditb_dim2 * diditb_dim1] = diditb[nd *
2624 diditb_dim2 * diditb_dim1] - bid1 * sou1 * sov1 -
2625 bid2 * sou2 * sov1 - bid3 * sou1 * sov2 - bid4 *
2626 sou2 * sov2;
2627 }
2628
2629/* L300: */
2630 }
2631/* L200: */
2632 }
2633/* L100: */
2634 }
2635 goto L9999;
2636
2637/* ------------------------------ The End -------------------------------
2638*/
2639
2640L9999:
2641 if (ibb >= 3) {
2642 AdvApp2Var_SysBase::mgsomsg_("MMA2CD1", 7L);
2643 }
2644 return 0;
2645} /* mma2cd1_ */
2646
2647//=======================================================================
2648//function : mma2cd2_
2649//purpose :
2650//=======================================================================
2651int mma2cd2_(integer *ndimen,
2652 integer *nbpntu,
2653 integer *nbpntv,
2654 doublereal *vrootl,
2655 integer *iordrv,
2656 doublereal *sotbv1,
2657 doublereal *sotbv2,
2658 doublereal *ditbv1,
2659 doublereal *ditbv2,
2660 doublereal *fpntab,
2661 doublereal *vhermt,
2662 doublereal *sosotb,
2663 doublereal *soditb,
2664 doublereal *disotb,
2665 doublereal *diditb)
2666
2667{
1ef32e96 2668 integer c__1 = 1;
7fd59977 2669 /* System generated locals */
2670 integer sotbv1_dim1, sotbv1_dim2, sotbv1_offset, sotbv2_dim1, sotbv2_dim2,
2671 sotbv2_offset, ditbv1_dim1, ditbv1_dim2, ditbv1_offset,
2672 ditbv2_dim1, ditbv2_dim2, ditbv2_offset, fpntab_dim1,
2673 fpntab_offset, vhermt_dim1, vhermt_offset, sosotb_dim1,
2674 sosotb_dim2, sosotb_offset, diditb_dim1, diditb_dim2,
2675 diditb_offset, soditb_dim1, soditb_dim2, soditb_offset,
2676 disotb_dim1, disotb_dim2, disotb_offset, i__1, i__2, i__3, i__4;
41194117 2677
7fd59977 2678 /* Local variables */
1ef32e96
RL
2679 integer ncfhv, nuroo, nvroo, ii, nd, jj, kk, ibb, jjm, jjp;
2680 doublereal bid1, bid2, bid3, bid4;
7fd59977 2681
2682/* **********************************************************************
2683*/
0d969553 2684/* FUNCTION : */
7fd59977 2685/* ---------- */
0d969553
Y
2686/* Discretisation on the parameters of polynoms of interpolation */
2687/* of constraints on 2 borders iso-V of order IORDRV. */
7fd59977 2688
0d969553
Y
2689
2690/* KEYWORDS : */
7fd59977 2691/* ----------- */
2692/* TOUS, AB_SPECIFI::CONTRAINTE&, DISCRETISATION, &POINT */
2693
7fd59977 2694
0d969553
Y
2695
2696/* INPUT ARGUMENTS : */
2697/* ------------------ */
2698/* NDIMEN: Dimension of the space. */
2699/* NBPNTU: Nb of INTERNAL parameters of discretisation by U. */
2700/* This is also the nb of root of Legendre polynom where discretization is done. */
2701/* UROOTL: Table of parameters of discretisation ON (-1,1) by U.
2702*/
2703/* NBPNTV: Nb of INTERNAL parameters of discretisation by V. */
2704/* This is also the nb of root of Legendre polynom where discretization is done. */
258ff83b 2705/* VROOTL: Table of discretization parameters on (-1,1) by V. */
0d969553
Y
2706/* IORDRV: Order of constraint imposed at the extremities of iso-V */
2707/* = 0, calculate the extremities of iso-V */
2708/* = 1, calculate, additionally, the 1st derivative in the direction of iso-V */
2709/* = 2, calculate, additionally, the 2nd derivative in the direction of iso-V */
2710/* SOTBV1: Table of NBPNTV/2 sums of 2 index points */
2711/* NBPNTV-II+1 and II, for II = 1, NBPNTV/2 on iso-V0. */
2712/* SOTBV2: Table of NBPNTV/2 sums of 2 index points */
2713/* NBPNTV-II+1 and II, for II = 1, NBPNTV/2 on iso-V1. */
2714/* DITBV1: Table of NBPNTV/2 differences of 2 index points */
2715/* NBPNTV-II+1 and II, for II = 1, NBPNTV/2 on iso-V0. */
2716/* DITBV2: Table of NBPNTV/2 differences of 2 index points */
2717/* NBPNTV-II+1 and II, for II = 1, NBPNTV/2 on iso-V1. */
2718/* SOSOTB: Preinitialized table (input/output argument). */
2719/* DISOTB: Preinitialized table (input/output argument). */
2720/* SODITB: Preinitialized table (input/output argument). */
2721/* DIDITB: Preinitialized table (input/output argument) */
2722
2723/* OUTPUT ARGUMENTS : */
7fd59977 2724/* ------------------- */
0d969553
Y
2725/* FPNTAB: Auxiliary table. */
2726/* VHERMT: Table of 2*(IORDRV+1) coeff. of 2*(IORDRV+1) polynoms of Hermite. */
2727/* SOSOTB: Table where the terms of constraints are added */
7fd59977 2728/* C(ui,vj) + C(ui,-vj) + C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
2729/* with ui and vj positive roots of the Legendre polynom */
2730/* of degree NBPNTU and NBPNTV respectively. */
2731/* DISOTB: Table where the terms of constraints are added */
7fd59977 2732/* C(ui,vj) + C(ui,-vj) - C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
2733/* with ui and vj positive roots of the polynom of Legendre */
2734/* of degree NBPNTU and NBPNTV respectively. */
2735/* SODITB: Table where the terms of constraints are added */
7fd59977 2736/* C(ui,vj) - C(ui,-vj) + C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
2737/* with ui and vj positive roots of the polynom of Legendre */
2738/* of degree NBPNTU and NBPNTV respectively. */
2739/* DIDITB: Table where the terms of constraints are added */
7fd59977 2740/* C(ui,vj) - C(ui,-vj) - C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
2741/* with ui and vj positive roots of the polynom of Legendre */
2742/* of degree NBPNTU and NBPNTV respectively. */
7fd59977 2743
0d969553 2744/* COMMONS USED : */
7fd59977 2745/* ---------------- */
2746
0d969553 2747/* REFERENCES CALLED : */
7fd59977 2748/* ----------------------- */
2749
0d969553 2750/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 2751/* ----------------------------------- */
2752
2753
7fd59977 2754/* > */
2755/* **********************************************************************
2756*/
2757
0d969553 2758/* Name of the routine */
7fd59977 2759
2760
2761 /* Parameter adjustments */
2762 diditb_dim1 = *nbpntu / 2 + 1;
2763 diditb_dim2 = *nbpntv / 2 + 1;
2764 diditb_offset = diditb_dim1 * diditb_dim2;
2765 diditb -= diditb_offset;
2766 disotb_dim1 = *nbpntu / 2;
2767 disotb_dim2 = *nbpntv / 2;
2768 disotb_offset = disotb_dim1 * (disotb_dim2 + 1) + 1;
2769 disotb -= disotb_offset;
2770 soditb_dim1 = *nbpntu / 2;
2771 soditb_dim2 = *nbpntv / 2;
2772 soditb_offset = soditb_dim1 * (soditb_dim2 + 1) + 1;
2773 soditb -= soditb_offset;
2774 sosotb_dim1 = *nbpntu / 2 + 1;
2775 sosotb_dim2 = *nbpntv / 2 + 1;
2776 sosotb_offset = sosotb_dim1 * sosotb_dim2;
2777 sosotb -= sosotb_offset;
2778 --vrootl;
2779 vhermt_dim1 = (*iordrv << 1) + 2;
2780 vhermt_offset = vhermt_dim1;
2781 vhermt -= vhermt_offset;
2782 fpntab_dim1 = *nbpntv;
2783 fpntab_offset = fpntab_dim1 + 1;
2784 fpntab -= fpntab_offset;
2785 ditbv2_dim1 = *nbpntu / 2 + 1;
2786 ditbv2_dim2 = *ndimen;
2787 ditbv2_offset = ditbv2_dim1 * (ditbv2_dim2 + 1);
2788 ditbv2 -= ditbv2_offset;
2789 ditbv1_dim1 = *nbpntu / 2 + 1;
2790 ditbv1_dim2 = *ndimen;
2791 ditbv1_offset = ditbv1_dim1 * (ditbv1_dim2 + 1);
2792 ditbv1 -= ditbv1_offset;
2793 sotbv2_dim1 = *nbpntu / 2 + 1;
2794 sotbv2_dim2 = *ndimen;
2795 sotbv2_offset = sotbv2_dim1 * (sotbv2_dim2 + 1);
2796 sotbv2 -= sotbv2_offset;
2797 sotbv1_dim1 = *nbpntu / 2 + 1;
2798 sotbv1_dim2 = *ndimen;
2799 sotbv1_offset = sotbv1_dim1 * (sotbv1_dim2 + 1);
2800 sotbv1 -= sotbv1_offset;
2801
2802 /* Function Body */
2803 ibb = AdvApp2Var_SysBase::mnfndeb_();
2804 if (ibb >= 3) {
2805 AdvApp2Var_SysBase::mgenmsg_("MMA2CD2", 7L);
2806 }
2807
0d969553 2808/* ------------------- Discretization of Hermit polynoms -----------
7fd59977 2809*/
2810
2811 ncfhv = (*iordrv + 1) << 1;
2812 i__1 = ncfhv;
2813 for (ii = 1; ii <= i__1; ++ii) {
2814 i__2 = *nbpntv;
2815 for (jj = 1; jj <= i__2; ++jj) {
2816 AdvApp2Var_MathBase::mmmpocur_(&ncfhv, &c__1, &ncfhv, &vhermt[ii * vhermt_dim1], &
2817 vrootl[jj], &fpntab[jj + ii * fpntab_dim1]);
2818/* L60: */
2819 }
2820/* L50: */
2821 }
2822
0d969553 2823/* ---- The discretizations of polynoms of constraints are subtracted ----
7fd59977 2824*/
2825
2826 nuroo = *nbpntu / 2;
2827 nvroo = *nbpntv / 2;
2828
2829 i__1 = *ndimen;
2830 for (nd = 1; nd <= i__1; ++nd) {
2831 i__2 = *iordrv + 1;
2832 for (ii = 1; ii <= i__2; ++ii) {
2833
2834 i__3 = nuroo;
2835 for (kk = 1; kk <= i__3; ++kk) {
2836 bid1 = sotbv1[kk + (nd + ii * sotbv1_dim2) * sotbv1_dim1];
2837 bid2 = sotbv2[kk + (nd + ii * sotbv2_dim2) * sotbv2_dim1];
2838 bid3 = ditbv1[kk + (nd + ii * ditbv1_dim2) * ditbv1_dim1];
2839 bid4 = ditbv2[kk + (nd + ii * ditbv2_dim2) * ditbv2_dim1];
2840 i__4 = nvroo;
2841 for (jj = 1; jj <= i__4; ++jj) {
2842 jjp = (*nbpntv + 1) / 2 + jj;
2843 jjm = nvroo - jj + 1;
2844 sosotb[kk + (jj + nd * sosotb_dim2) * sosotb_dim1] =
2845 sosotb[kk + (jj + nd * sosotb_dim2) * sosotb_dim1]
2846 - bid1 * (fpntab[jjp + ((ii << 1) - 1) *
2847 fpntab_dim1] + fpntab[jjm + ((ii << 1) - 1) *
2848 fpntab_dim1]) - bid2 * (fpntab[jjp + (ii << 1) *
2849 fpntab_dim1] + fpntab[jjm + (ii << 1) *
2850 fpntab_dim1]);
2851 disotb[kk + (jj + nd * disotb_dim2) * disotb_dim1] =
2852 disotb[kk + (jj + nd * disotb_dim2) * disotb_dim1]
2853 - bid3 * (fpntab[jjp + ((ii << 1) - 1) *
2854 fpntab_dim1] + fpntab[jjm + ((ii << 1) - 1) *
2855 fpntab_dim1]) - bid4 * (fpntab[jjp + (ii << 1) *
2856 fpntab_dim1] + fpntab[jjm + (ii << 1) *
2857 fpntab_dim1]);
2858 soditb[kk + (jj + nd * soditb_dim2) * soditb_dim1] =
2859 soditb[kk + (jj + nd * soditb_dim2) * soditb_dim1]
2860 - bid1 * (fpntab[jjp + ((ii << 1) - 1) *
2861 fpntab_dim1] - fpntab[jjm + ((ii << 1) - 1) *
2862 fpntab_dim1]) - bid2 * (fpntab[jjp + (ii << 1) *
2863 fpntab_dim1] - fpntab[jjm + (ii << 1) *
2864 fpntab_dim1]);
2865 diditb[kk + (jj + nd * diditb_dim2) * diditb_dim1] =
2866 diditb[kk + (jj + nd * diditb_dim2) * diditb_dim1]
2867 - bid3 * (fpntab[jjp + ((ii << 1) - 1) *
2868 fpntab_dim1] - fpntab[jjm + ((ii << 1) - 1) *
2869 fpntab_dim1]) - bid4 * (fpntab[jjp + (ii << 1) *
2870 fpntab_dim1] - fpntab[jjm + (ii << 1) *
2871 fpntab_dim1]);
2872/* L400: */
2873 }
2874/* L300: */
2875 }
2876/* L200: */
2877 }
2878
0d969553
Y
2879/* ------------ Case when the discretization is done only on the roots */
2880/* ---------- of Legendre polynom of uneven degree, 0 is root */
2881
7fd59977 2882
2883 if (*nbpntv % 2 == 1) {
2884 i__2 = *iordrv + 1;
2885 for (ii = 1; ii <= i__2; ++ii) {
2886 i__3 = nuroo;
2887 for (kk = 1; kk <= i__3; ++kk) {
2888 bid1 = sotbv1[kk + (nd + ii * sotbv1_dim2) * sotbv1_dim1]
2889 * fpntab[nvroo + 1 + ((ii << 1) - 1) *
2890 fpntab_dim1] + sotbv2[kk + (nd + ii * sotbv2_dim2)
2891 * sotbv2_dim1] * fpntab[nvroo + 1 + (ii << 1) *
2892 fpntab_dim1];
2893 sosotb[kk + nd * sosotb_dim2 * sosotb_dim1] -= bid1;
2894 bid2 = ditbv1[kk + (nd + ii * ditbv1_dim2) * ditbv1_dim1]
2895 * fpntab[nvroo + 1 + ((ii << 1) - 1) *
2896 fpntab_dim1] + ditbv2[kk + (nd + ii * ditbv2_dim2)
2897 * ditbv2_dim1] * fpntab[nvroo + 1 + (ii << 1) *
2898 fpntab_dim1];
2899 diditb[kk + nd * diditb_dim2 * diditb_dim1] -= bid2;
2900/* L550: */
2901 }
2902/* L500: */
2903 }
2904 }
2905
2906 if (*nbpntu % 2 == 1) {
2907 i__2 = *iordrv + 1;
2908 for (ii = 1; ii <= i__2; ++ii) {
2909 i__3 = nvroo;
2910 for (jj = 1; jj <= i__3; ++jj) {
2911 jjp = (*nbpntv + 1) / 2 + jj;
2912 jjm = nvroo - jj + 1;
2913 bid1 = sotbv1[(nd + ii * sotbv1_dim2) * sotbv1_dim1] * (
2914 fpntab[jjp + ((ii << 1) - 1) * fpntab_dim1] +
2915 fpntab[jjm + ((ii << 1) - 1) * fpntab_dim1]) +
2916 sotbv2[(nd + ii * sotbv2_dim2) * sotbv2_dim1] * (
2917 fpntab[jjp + (ii << 1) * fpntab_dim1] + fpntab[
2918 jjm + (ii << 1) * fpntab_dim1]);
2919 sosotb[(jj + nd * sosotb_dim2) * sosotb_dim1] -= bid1;
2920 bid2 = sotbv1[(nd + ii * sotbv1_dim2) * sotbv1_dim1] * (
2921 fpntab[jjp + ((ii << 1) - 1) * fpntab_dim1] -
2922 fpntab[jjm + ((ii << 1) - 1) * fpntab_dim1]) +
2923 sotbv2[(nd + ii * sotbv2_dim2) * sotbv2_dim1] * (
2924 fpntab[jjp + (ii << 1) * fpntab_dim1] - fpntab[
2925 jjm + (ii << 1) * fpntab_dim1]);
2926 diditb[jj + nd * diditb_dim2 * diditb_dim1] -= bid2;
2927/* L650: */
2928 }
2929/* L600: */
2930 }
2931 }
2932
2933 if (*nbpntu % 2 == 1 && *nbpntv % 2 == 1) {
2934 i__2 = *iordrv + 1;
2935 for (ii = 1; ii <= i__2; ++ii) {
2936 bid1 = sotbv1[(nd + ii * sotbv1_dim2) * sotbv1_dim1] * fpntab[
2937 nvroo + 1 + ((ii << 1) - 1) * fpntab_dim1] + sotbv2[(
2938 nd + ii * sotbv2_dim2) * sotbv2_dim1] * fpntab[nvroo
2939 + 1 + (ii << 1) * fpntab_dim1];
2940 sosotb[nd * sosotb_dim2 * sosotb_dim1] -= bid1;
2941/* L700: */
2942 }
2943 }
2944
2945/* L100: */
2946 }
2947 goto L9999;
2948
2949/* ------------------------------ The End -------------------------------
2950*/
2951
2952L9999:
2953 if (ibb >= 3) {
2954 AdvApp2Var_SysBase::mgsomsg_("MMA2CD2", 7L);
2955 }
2956 return 0;
2957} /* mma2cd2_ */
2958
2959//=======================================================================
2960//function : mma2cd3_
2961//purpose :
2962//=======================================================================
2963int mma2cd3_(integer *ndimen,
2964 integer *nbpntu,
2965 doublereal *urootl,
2966 integer *nbpntv,
2967 integer *iordru,
2968 doublereal *sotbu1,
2969 doublereal *sotbu2,
2970 doublereal *ditbu1,
2971 doublereal *ditbu2,
2972 doublereal *fpntab,
2973 doublereal *uhermt,
2974 doublereal *sosotb,
2975 doublereal *soditb,
2976 doublereal *disotb,
2977 doublereal *diditb)
2978
2979{
1ef32e96 2980 integer c__1 = 1;
41194117 2981
7fd59977 2982 /* System generated locals */
2983 integer sotbu1_dim1, sotbu1_dim2, sotbu1_offset, sotbu2_dim1, sotbu2_dim2,
2984 sotbu2_offset, ditbu1_dim1, ditbu1_dim2, ditbu1_offset,
2985 ditbu2_dim1, ditbu2_dim2, ditbu2_offset, fpntab_dim1,
2986 fpntab_offset, uhermt_dim1, uhermt_offset, sosotb_dim1,
2987 sosotb_dim2, sosotb_offset, diditb_dim1, diditb_dim2,
2988 diditb_offset, soditb_dim1, soditb_dim2, soditb_offset,
2989 disotb_dim1, disotb_dim2, disotb_offset, i__1, i__2, i__3, i__4;
2990
2991 /* Local variables */
1ef32e96
RL
2992 integer ncfhu, nuroo, nvroo, ii, nd, jj, kk, ibb, kkm, kkp;
2993 doublereal bid1, bid2, bid3, bid4;
7fd59977 2994
2995/* **********************************************************************
2996*/
0d969553 2997/* FUNCTION : */
7fd59977 2998/* ---------- */
0d969553
Y
2999/* Discretisation on the parameters of polynoms of interpolation */
3000/* of constraints on 2 borders iso-U of order IORDRU. */
7fd59977 3001
0d969553
Y
3002
3003/* KEYWORDS : */
7fd59977 3004/* ----------- */
3005/* TOUS, AB_SPECIFI::CONTRAINTE&, DISCRETISATION, &POINT */
3006
0d969553 3007/* INPUT ARGUMENTS : */
7fd59977 3008/* ------------------ */
0d969553
Y
3009/* NDIMEN: Dimension of the space. */
3010/* NBPNTU: Nb of INTERNAL parameters of discretisation by U. */
3011/* This is also the nb of root of Legendre polynom where discretization is done. */
3012/* UROOTL: Table of parameters of discretisation ON (-1,1) by U.
3013*/
3014/* NBPNTV: Nb of INTERNAL parameters of discretisation by V. */
3015/* This is also the nb of root of Legendre polynom where discretization is done. */
3016/* IORDRV: Order of constraint imposed at the extremities of iso-V */
3017/* = 0, calculate the extremities of iso-V */
3018/* = 1, calculate, additionally, the 1st derivative in the direction of iso-V */
3019/* = 2, calculate, additionally, the 2nd derivative in the direction of iso-V */
3020/* SOTBU1: Table of NBPNTU/2 sums of 2 index points */
3021/* NBPNTU-II+1 and II, for II = 1, NBPNTU/2 on iso-V0. */
3022/* SOTBU2: Table of NBPNTV/2 sums of 2 index points */
3023/* NBPNTU-II+1 and II, for II = 1, NBPNTU/2 on iso-V1. */
3024/* DITBU1: Table of NBPNTU/2 differences of 2 index points */
3025/* NBPNTU-II+1 and II, for II = 1, NBPNTU/2 on iso-V0. */
3026/* DITBU2: Table of NBPNTU/2 differences of 2 index points */
3027/* NBPNTU-II+1 and II, for II = 1, NBPNTU/2 on iso-V1. */
3028/* SOSOTB: Preinitialized table (input/output argument). */
3029/* DISOTB: Preinitialized table (input/output argument). */
3030/* SODITB: Preinitialized table (input/output argument). */
3031/* DIDITB: Preinitialized table (input/output argument) */
3032
3033/* OUTPUT ARGUMENTS : */
7fd59977 3034/* ------------------- */
0d969553
Y
3035/* FPNTAB: Auxiliary table. */
3036/* UHERMT: Table of 2*(IORDRU+1) coeff. of 2*(IORDRU+1) polynoms of Hermite. */
3037/* SOSOTB: Table where the terms of constraints are added */
7fd59977 3038/* C(ui,vj) + C(ui,-vj) + C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
3039/* with ui and vj positive roots of the Legendre polynom */
3040/* of degree NBPNTU and NBPNTV respectively. */
3041/* DISOTB: Table where the terms of constraints are added */
7fd59977 3042/* C(ui,vj) + C(ui,-vj) - C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
3043/* with ui and vj positive roots of the polynom of Legendre */
3044/* of degree NBPNTU and NBPNTV respectively. */
3045/* SODITB: Table where the terms of constraints are added */
7fd59977 3046/* C(ui,vj) - C(ui,-vj) + C(-ui,vj) - C(-ui,-vj) */
0d969553
Y
3047/* with ui and vj positive roots of the polynom of Legendre */
3048/* of degree NBPNTU and NBPNTV respectively. */
3049/* DIDITB: Table where the terms of constraints are added */
7fd59977 3050/* C(ui,vj) - C(ui,-vj) - C(-ui,vj) + C(-ui,-vj) */
0d969553
Y
3051/* with ui and vj positive roots of the polynom of Legendre */
3052/* of degree NBPNTU and NBPNTV respectively. */
7fd59977 3053
0d969553 3054/* COMMONS USED : */
7fd59977 3055/* ---------------- */
3056
0d969553 3057/* REFERENCES CALLED : */
7fd59977 3058/* ----------------------- */
3059
0d969553 3060/* DESCRIPTION/NOTES/LIMITATIONS : */
7fd59977 3061/* ----------------------------------- */
3062
7fd59977 3063/* $ HISTORIQUE DES MODIFICATIONS : */
3064/* -------------------------------- */
3065/* 08-08-1991: RBD; Creation. */
3066/* > */
3067/* **********************************************************************
3068*/
3069
0d969553 3070/* Name of the routine */
7fd59977 3071
3072
3073 /* Parameter adjustments */
3074 --urootl;
3075 diditb_dim1 = *nbpntu / 2 + 1;
3076 diditb_dim2 = *nbpntv / 2 + 1;
3077 diditb_offset = diditb_dim1 * diditb_dim2;
3078 diditb -= diditb_offset;
3079 disotb_dim1 = *nbpntu / 2;
3080 disotb_dim2 = *nbpntv / 2;
3081 disotb_offset = disotb_dim1 * (disotb_dim2 + 1) + 1;
3082 disotb -= disotb_offset;
3083 soditb_dim1 = *nbpntu / 2;
3084 soditb_dim2 = *nbpntv / 2;
3085 soditb_offset = soditb_dim1 * (soditb_dim2 + 1) + 1;
3086 soditb -= soditb_offset;
3087 sosotb_dim1 = *nbpntu / 2 + 1;
3088 sosotb_dim2 = *nbpntv / 2 + 1;
3089 sosotb_offset = sosotb_dim1 * sosotb_dim2;
3090 sosotb -= sosotb_offset;
3091 uhermt_dim1 = (*iordru << 1) + 2;
3092 uhermt_offset = uhermt_dim1;
3093 uhermt -= uhermt_offset;
3094 fpntab_dim1 = *nbpntu;
3095 fpntab_offset = fpntab_dim1 + 1;
3096 fpntab -= fpntab_offset;
3097 ditbu2_dim1 = *nbpntv / 2 + 1;
3098 ditbu2_dim2 = *ndimen;
3099 ditbu2_offset = ditbu2_dim1 * (ditbu2_dim2 + 1);
3100 ditbu2 -= ditbu2_offset;
3101 ditbu1_dim1 = *nbpntv / 2 + 1;
3102 ditbu1_dim2 = *ndimen;
3103 ditbu1_offset = ditbu1_dim1 * (ditbu1_dim2 + 1);
3104 ditbu1 -= ditbu1_offset;
3105 sotbu2_dim1 = *nbpntv / 2 + 1;
3106 sotbu2_dim2 = *ndimen;
3107 sotbu2_offset = sotbu2_dim1 * (sotbu2_dim2 + 1);
3108 sotbu2 -= sotbu2_offset;
3109 sotbu1_dim1 = *nbpntv / 2 + 1;
3110 sotbu1_dim2 = *ndimen;
3111 sotbu1_offset = sotbu1_dim1 * (sotbu1_dim2 + 1);
3112 sotbu1 -= sotbu1_offset;
3113
3114 /* Function Body */
3115 ibb = AdvApp2Var_SysBase::mnfndeb_();
3116 if (ibb >= 3) {
3117 AdvApp2Var_SysBase::mgenmsg_("MMA2CD3", 7L);
3118 }
3119
0d969553 3120/* ------------------- Discretization of polynoms of Hermit -----------
7fd59977 3121*/
3122
3123 ncfhu = (*iordru + 1) << 1;
3124 i__1 = ncfhu;
3125 for (ii = 1; ii <= i__1; ++ii) {
3126 i__2 = *nbpntu;
3127 for (kk = 1; kk <= i__2; ++kk) {
3128 AdvApp2Var_MathBase::mmmpocur_(&ncfhu,
3129 &c__1,
3130 &ncfhu,
3131 &uhermt[ii * uhermt_dim1],
3132 &urootl[kk],
3133 &fpntab[kk + ii * fpntab_dim1]);
3134/* L60: */
3135 }
3136/* L50: */
3137 }
3138
0d969553 3139/* ---- The discretizations of polynoms of constraints are subtracted ----
7fd59977 3140*/
3141
3142 nvroo = *nbpntv / 2;
3143 nuroo = *nbpntu / 2;
3144
3145 i__1 = *ndimen;
3146 for (nd = 1; nd <= i__1; ++nd) {
3147 i__2 = *iordru + 1;
3148 for (ii = 1; ii <= i__2; ++ii) {
3149
3150 i__3 = nvroo;
3151 for (jj = 1; jj <= i__3; ++jj) {
3152 bid1 = sotbu1[jj + (nd + ii * sotbu1_dim2) * sotbu1_dim1];
3153 bid2 = sotbu2[jj + (nd + ii * sotbu2_dim2) * sotbu2_dim1];
3154 bid3 = ditbu1[jj + (nd + ii * ditbu1_dim2) * ditbu1_dim1];
3155 bid4 = ditbu2[jj + (nd + ii * ditbu2_dim2) * ditbu2_dim1];
3156 i__4 = nuroo;
3157 for (kk = 1; kk <= i__4; ++kk) {
3158 kkp = (*nbpntu + 1) / 2 + kk;
3159 kkm = nuroo - kk + 1;
3160 sosotb[kk + (jj + nd * sosotb_dim2) * sosotb_dim1] =
3161 sosotb[kk + (jj + nd * sosotb_dim2) * sosotb_dim1]
3162 - bid1 * (fpntab[kkp + ((ii << 1) - 1) *
3163 fpntab_dim1] + fpntab[kkm + ((ii << 1) - 1) *
3164 fpntab_dim1]) - bid2 * (fpntab[kkp + (ii << 1) *
3165 fpntab_dim1] + fpntab[kkm + (ii << 1) *
3166 fpntab_dim1]);
3167 disotb[kk + (jj + nd * disotb_dim2) * disotb_dim1] =
3168 disotb[kk + (jj + nd * disotb_dim2) * disotb_dim1]
3169 - bid1 * (fpntab[kkp + ((ii << 1) - 1) *
3170 fpntab_dim1] - fpntab[kkm + ((ii << 1) - 1) *
3171 fpntab_dim1]) - bid2 * (fpntab[kkp + (ii << 1) *
3172 fpntab_dim1] - fpntab[kkm + (ii << 1) *
3173 fpntab_dim1]);
3174 soditb[kk + (jj + nd * soditb_dim2) * soditb_dim1] =
3175 soditb[kk + (jj + nd * soditb_dim2) * soditb_dim1]
3176 - bid3 * (fpntab[kkp + ((ii << 1) - 1) *
3177 fpntab_dim1] + fpntab[kkm + ((ii << 1) - 1) *
3178 fpntab_dim1]) - bid4 * (fpntab[kkp + (ii << 1) *
3179 fpntab_dim1] + fpntab[kkm + (ii << 1) *
3180 fpntab_dim1]);
3181 diditb[kk + (jj + nd * diditb_dim2) * diditb_dim1] =
3182 diditb[kk + (jj + nd * diditb_dim2) * diditb_dim1]
3183 - bid3 * (fpntab[kkp + ((ii << 1) - 1) *
3184 fpntab_dim1] - fpntab[kkm + ((ii << 1) - 1) *
3185 fpntab_dim1]) - bid4 * (fpntab[kkp + (ii << 1) *
3186 fpntab_dim1] - fpntab[kkm + (ii << 1) *
3187 fpntab_dim1]);
3188/* L400: */
3189 }
3190/* L300: */
3191 }
3192/* L200: */
3193 }
3194
0d969553
Y
3195/* ------------ Case when the discretization is done only on the roots */
3196/* ---------- of Legendre polynom of uneven degree, 0 is root */
3197
3198
7fd59977 3199
3200 if (*nbpntu % 2 == 1) {
3201 i__2 = *iordru + 1;
3202 for (ii = 1; ii <= i__2; ++ii) {
3203 i__3 = nvroo;
3204 for (jj = 1; jj <= i__3; ++jj) {
3205 bid1 = sotbu1[jj + (nd + ii * sotbu1_dim2) * sotbu1_dim1]
3206 * fpntab[nuroo + 1 + ((ii << 1) - 1) *
3207 fpntab_dim1] + sotbu2[jj + (nd + ii * sotbu2_dim2)
3208 * sotbu2_dim1] * fpntab[nuroo + 1 + (ii << 1) *
3209 fpntab_dim1];
3210 sosotb[(jj + nd * sosotb_dim2) * sosotb_dim1] -= bid1;
3211 bid2 = ditbu1[jj + (nd + ii * ditbu1_dim2) * ditbu1_dim1]
3212 * fpntab[nuroo + 1 + ((ii << 1) - 1) *
3213 fpntab_dim1] + ditbu2[jj + (nd + ii * ditbu2_dim2)
3214 * ditbu2_dim1] * fpntab[nuroo + 1 + (ii << 1) *
3215 fpntab_dim1];
3216 diditb[(jj + nd * diditb_dim2) * diditb_dim1] -= bid2;
3217/* L550: */
3218 }
3219/* L500: */
3220 }
3221 }
3222
3223 if (*nbpntv % 2 == 1) {
3224 i__2 = *iordru + 1;
3225 for (ii = 1; ii <= i__2; ++ii) {
3226 i__3 = nuroo;
3227 for (kk = 1; kk <= i__3; ++kk) {
3228 kkp = (*nbpntu + 1) / 2 + kk;
3229 kkm = nuroo - kk + 1;
3230 bid1 = sotbu1[(nd + ii * sotbu1_dim2) * sotbu1_dim1] * (
3231 fpntab[kkp + ((ii << 1) - 1) * fpntab_dim1] +
3232 fpntab[kkm + ((ii << 1) - 1) * fpntab_dim1]) +
3233 sotbu2[(nd + ii * sotbu2_dim2) * sotbu2_dim1] * (
3234 fpntab[kkp + (ii << 1) * fpntab_dim1] + fpntab[
3235 kkm + (ii << 1) * fpntab_dim1]);
3236 sosotb[kk + nd * sosotb_dim2 * sosotb_dim1] -= bid1;
3237 bid2 = sotbu1[(nd + ii * sotbu1_dim2) * sotbu1_dim1] * (
3238 fpntab[kkp + ((ii << 1) - 1) * fpntab_dim1] -
3239 fpntab[kkm + ((ii << 1) - 1) * fpntab_dim1]) +
3240 sotbu2[(nd + ii * sotbu2_dim2) * sotbu2_dim1] * (
3241 fpntab[kkp + (ii << 1) * fpntab_dim1] - fpntab[
3242 kkm + (ii << 1) * fpntab_dim1]);
3243 diditb[kk + nd * diditb_dim2 * diditb_dim1] -= bid2;
3244/* L650: */
3245 }
3246/* L600: */
3247 }
3248 }
3249
3250 if (*nbpntu % 2 == 1 && *nbpntv % 2 == 1) {
3251 i__2 = *iordru + 1;
3252 for (ii = 1; ii <= i__2; ++ii) {
3253 bid1 = sotbu1[(nd + ii * sotbu1_dim2) * sotbu1_dim1] * fpntab[
3254 nuroo + 1 + ((ii << 1) - 1) * fpntab_dim1] + sotbu2[(
3255 nd + ii * sotbu2_dim2) * sotbu2_dim1] * fpntab[nuroo
3256 + 1 + (ii << 1) * fpntab_dim1];
3257 sosotb[nd * sosotb_dim2 * sosotb_dim1] -= bid1;
3258/* L700: */
3259 }
3260 }
3261
3262/* L100: */
3263 }
3264 goto L9999;
3265
3266/* ------------------------------ The End -------------------------------
3267*/
3268
3269L9999:
3270 if (ibb >= 3) {
3271 AdvApp2Var_SysBase::mgsomsg_("MMA2CD3", 7L);
3272 }
3273 return 0;
3274} /* mma2cd3_ */
3275
3276//=======================================================================
3277//function : mma2cdi_
3278//purpose :
3279//=======================================================================
3280int AdvApp2Var_ApproxF2var::mma2cdi_( integer *ndimen,
3281 integer *nbpntu,
3282 doublereal *urootl,
3283 integer *nbpntv,
3284 doublereal *vrootl,
3285 integer *iordru,
3286 integer *iordrv,
3287 doublereal *contr1,
3288 doublereal *contr2,
3289 doublereal *contr3,
3290 doublereal *contr4,
3291 doublereal *sotbu1,
3292 doublereal *sotbu2,
3293 doublereal *ditbu1,
3294 doublereal *ditbu2,
3295 doublereal *sotbv1,
3296 doublereal *sotbv2,
3297 doublereal *ditbv1,
3298 doublereal *ditbv2,
3299 doublereal *sosotb,
3300 doublereal *soditb,
3301 doublereal *disotb,
3302 doublereal *diditb,
3303 integer *iercod)
3304
3305{
1ef32e96 3306 integer c__8 = 8;
7fd59977 3307
3308 /* System generated locals */
3309 integer contr1_dim1, contr1_dim2, contr1_offset, contr2_dim1, contr2_dim2,
3310 contr2_offset, contr3_dim1, contr3_dim2, contr3_offset,
3311 contr4_dim1, contr4_dim2, contr4_offset, sosotb_dim1, sosotb_dim2,
3312 sosotb_offset, diditb_dim1, diditb_dim2, diditb_offset,
3313 soditb_dim1, soditb_dim2, soditb_offset, disotb_dim1, disotb_dim2,