454 {
455 int j, j1, k, k1, ke, iloc, err = 0;
456 float dt1, xl0, xl1, xl2, x, y, qcx, qcy, airne;
457 creal gconvectke;
458
459 if (!__quadra)
460 {
461 if (i1 == ns)
462 {
463 creal aux = gconvect[i1];
464
465 delete [] gconvect;gconvect = NULL;
466 return aux;
467 }
468 else if (i1 > 0)
469 return gconvect[i1];
470 else
471 {
472 gconvect = new creal[ns];
473 rhsQuadra = 1;
474 for (ke = 0; ke < ne; ke++)
475 {
476 airne = 0.F;
477 if (k1 = Tleft[ke], k1 >= 0)
478 airne += area[k1];
479 if (k = Tright[ke], k >= 0)
480 airne += area[k];
481 else
482 k = k1;
483 qcx = (q[me[k][0]][0] + q[me[k][1]][0] + q[me[k][2]][0]) / 3;
484 qcy = (q[me[k][0]][1] + q[me[k][1]][1] + q[me[k][2]][1]) / 3;
485 j = lowV[ke];
486 j1 = highV[ke];
487 x = qcx + 0.999F * ((q[j][0] + q[j1][0]) / 2 - qcx);
488 y = qcy + 0.999F * ((q[j][1] + q[j1][1]) / 2 - qcy);
489 dt1 = dt;
490 k1 = k;
491 if (0 < xtoX (u1, u2, &dt1, &x, &y, &k1))
492 err++;
493 if (barycoor (x, y, k1, &xl0, &xl1, &xl2))
494 err += 3;
495 gconvectke = f[me[k1][0]] * xl0 + f[me[k1][1]] * xl1 + f[me[k1][2]] * xl2;
496 gconvect[lowV[ke]] += airne * gconvectke / 6.F;
497 gconvect[highV[ke]] += airne * gconvectke / 6.F;
498 }
499 return gconvect[0];
500 }
501 }
502 else
503
504 {
505 if (i1 == 3 * nt - 1)
506 {
507 creal aux = gconvect[i1];
508
509 delete [] gconvect;gconvect = NULL;
510 return aux;
511 }
512 else if (i1 > 0)
513 return gconvect[i1];
514 else
515
516 {
517 creal gloc[3];
518
519 gconvect = new creal[ns];
520 for (k = 0; k < nt; k++)
521 {
522 qcx = (q[me[k][0]][0] + q[me[k][1]][0] + q[me[k][2]][0]) / 3;
523 qcy = (q[me[k][0]][1] + q[me[k][1]][1] + q[me[k][2]][1]) / 3;
524 for (iloc = 0; iloc < 3; iloc++)
525 {
526 j = me[k][iloc];
527 j1 = me[k][next[iloc]];
528 x = qcx + 0.999F * ((q[j][0] + q[j1][0]) / 2 - qcx);
529 y = qcy + 0.999F * ((q[j][1] + q[j1][1]) / 2 - qcy);
530 dt1 = dt;
531 k1 = k;
532 if (0 < xtoX (u1, u2, &dt1, &x, &y, &k1))
533 err++;
534 if (barycoor (x, y, k1, &xl0, &xl1, &xl2))
535 err += 3;
536 gloc[next[iloc + 1]] = f[3 * k1] * xl0 + f[3 * k1 + 1] * xl1 + f[3 * k1 + 2] * xl2;
537 }
538 for (iloc = 0; iloc < 3; iloc++)
539 gconvect[3 * k + iloc] = gloc[next[iloc]] + gloc[next[iloc + 1]] - gloc[iloc];
540 }
541 return gconvect[0];
542 }
543 }
544}