| [70c4392] | 1 | /*
|
|---|
| 2 | A simple 2D hydro code
|
|---|
| 3 | (C) Romain Teyssier : CEA/IRFU -- original F90 code
|
|---|
| 4 | (C) Pierre-Francois Lavallee : IDRIS -- original F90 code
|
|---|
| 5 | (C) Guillaume Colin de Verdiere : CEA/DAM -- for the C version
|
|---|
| 6 | */
|
|---|
| 7 | /*
|
|---|
| 8 |
|
|---|
| 9 | This software is governed by the CeCILL license under French law and
|
|---|
| 10 | abiding by the rules of distribution of free software. You can use,
|
|---|
| 11 | modify and/ or redistribute the software under the terms of the CeCILL
|
|---|
| 12 | license as circulated by CEA, CNRS and INRIA at the following URL
|
|---|
| 13 | "http://www.cecill.info".
|
|---|
| 14 |
|
|---|
| 15 | As a counterpart to the access to the source code and rights to copy,
|
|---|
| 16 | modify and redistribute granted by the license, users are provided only
|
|---|
| 17 | with a limited warranty and the software's author, the holder of the
|
|---|
| 18 | economic rights, and the successive licensors have only limited
|
|---|
| 19 | liability.
|
|---|
| 20 |
|
|---|
| 21 | In this respect, the user's attention is drawn to the risks associated
|
|---|
| 22 | with loading, using, modifying and/or developing or reproducing the
|
|---|
| 23 | software by the user in light of its specific status of free software,
|
|---|
| 24 | that may mean that it is complicated to manipulate, and that also
|
|---|
| 25 | therefore means that it is reserved for developers and experienced
|
|---|
| 26 | professionals having in-depth computer knowledge. Users are therefore
|
|---|
| 27 | encouraged to load and test the software's suitability as regards their
|
|---|
| 28 | requirements in conditions enabling the security of their systems and/or
|
|---|
| 29 | data to be ensured and, more generally, to use and operate it in the
|
|---|
| 30 | same conditions as regards security.
|
|---|
| 31 |
|
|---|
| 32 | The fact that you are presently reading this means that you have had
|
|---|
| 33 | knowledge of the CeCILL license and that you accept its terms.
|
|---|
| 34 |
|
|---|
| 35 | */
|
|---|
| 36 |
|
|---|
| 37 | #include <stdlib.h>
|
|---|
| 38 | #include <unistd.h>
|
|---|
| 39 | #include <math.h>
|
|---|
| 40 | #include <stdio.h>
|
|---|
| 41 |
|
|---|
| 42 | #include "parametres.h"
|
|---|
| 43 | #include "utils.h"
|
|---|
| 44 | #include "slope.h"
|
|---|
| 45 | #include "perfcnt.h"
|
|---|
| 46 |
|
|---|
| 47 | #ifndef HMPP
|
|---|
| 48 |
|
|---|
| 49 | #define DABS(x) (real_t) fabs((x))
|
|---|
| 50 |
|
|---|
| 51 | void
|
|---|
| 52 | slope(const int n,
|
|---|
| 53 | const int Hnvar,
|
|---|
| 54 | const int Hnxyt,
|
|---|
| 55 | const real_t Hslope_type,
|
|---|
| 56 | const int slices, const int Hstep, real_t q[Hnvar][Hstep][Hnxyt], real_t dq[Hnvar][Hstep][Hnxyt]) {
|
|---|
| 57 | int nbv, i, ijmin, ijmax, s;
|
|---|
| 58 | // long ihvwin, ihvwimn, ihvwipn;
|
|---|
| 59 | // #define IHVW(i, v) ((i) + (v) * Hnxyt)
|
|---|
| 60 |
|
|---|
| 61 | WHERE("slope");
|
|---|
| 62 | ijmin = 0;
|
|---|
| 63 | ijmax = n;
|
|---|
| 64 |
|
|---|
| 65 | // #define OLDSTYLE
|
|---|
| 66 |
|
|---|
| 67 | #pragma omp parallel for private(nbv, s, i) shared(dq) COLLAPSE
|
|---|
| 68 | for (s = 0; s < slices; s++) {
|
|---|
| 69 | for (nbv = 0; nbv < Hnvar; nbv++) {
|
|---|
| 70 | #pragma ivdep
|
|---|
| 71 | for (i = ijmin + 1; i < ijmax - 1; i++) {
|
|---|
| 72 | real_t dlft, drgt, dcen, dsgn, slop, dlim;
|
|---|
| 73 | int llftrgt = 0;
|
|---|
| 74 | real_t t1;
|
|---|
| 75 | dlft = Hslope_type * (q[nbv][s][i] - q[nbv][s][i - 1]);
|
|---|
| 76 | drgt = Hslope_type * (q[nbv][s][i + 1] - q[nbv][s][i]);
|
|---|
| 77 | dcen = half * (dlft + drgt) / Hslope_type;
|
|---|
| 78 | dsgn = (dcen > 0) ? (real_t) 1.0 : (real_t) -1.0; // sign(one, dcen);
|
|---|
| 79 | #ifdef OLDSTYLE
|
|---|
| 80 | slop = fmin(fabs(dlft), fabs(drgt));
|
|---|
| 81 | dlim = slop;
|
|---|
| 82 | if ((dlft * drgt) <= zero) {
|
|---|
| 83 | dlim = zero;
|
|---|
| 84 | }
|
|---|
| 85 | dq[nbv][s][i] = dsgn * fmin(dlim, fabs(dcen));
|
|---|
| 86 | #else
|
|---|
| 87 | llftrgt = ((dlft * drgt) <= zero);
|
|---|
| 88 | t1 = fmin(fabs(dlft), fabs(drgt));
|
|---|
| 89 | dq[nbv][s][i] = dsgn * fmin((1 - llftrgt) * t1, fabs(dcen));
|
|---|
| 90 | #endif
|
|---|
| 91 | }
|
|---|
| 92 | }
|
|---|
| 93 | }
|
|---|
| 94 | {
|
|---|
| 95 | int nops = Hnvar * slices * ((ijmax - 1) - (ijmin + 1));
|
|---|
| 96 | FLOPS(8 * nops, 1 * nops, 6 * nops, 0 * nops);
|
|---|
| 97 | }
|
|---|
| 98 | } // slope
|
|---|
| 99 |
|
|---|
| 100 | #undef IHVW
|
|---|
| 101 |
|
|---|
| 102 | #endif /* HMPP */
|
|---|
| 103 | //EOF
|
|---|