| 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 |
|
|---|
| 10 | This software is governed by the CeCILL license under French law and
|
|---|
| 11 | abiding by the rules of distribution of free software. You can use,
|
|---|
| 12 | modify and/ or redistribute the software under the terms of the CeCILL
|
|---|
| 13 | license as circulated by CEA, CNRS and INRIA at the following URL
|
|---|
| 14 | "http://www.cecill.info".
|
|---|
| 15 |
|
|---|
| 16 | As a counterpart to the access to the source code and rights to copy,
|
|---|
| 17 | modify and redistribute granted by the license, users are provided only
|
|---|
| 18 | with a limited warranty and the software's author, the holder of the
|
|---|
| 19 | economic rights, and the successive licensors have only limited
|
|---|
| 20 | liability.
|
|---|
| 21 |
|
|---|
| 22 | In this respect, the user's attention is drawn to the risks associated
|
|---|
| 23 | with loading, using, modifying and/or developing or reproducing the
|
|---|
| 24 | software by the user in light of its specific status of free software,
|
|---|
| 25 | that may mean that it is complicated to manipulate, and that also
|
|---|
| 26 | therefore means that it is reserved for developers and experienced
|
|---|
| 27 | professionals having in-depth computer knowledge. Users are therefore
|
|---|
| 28 | encouraged to load and test the software's suitability as regards their
|
|---|
| 29 | requirements in conditions enabling the security of their systems and/or
|
|---|
| 30 | data to be ensured and, more generally, to use and operate it in the
|
|---|
| 31 | same conditions as regards security.
|
|---|
| 32 |
|
|---|
| 33 | The fact that you are presently reading this means that you have had
|
|---|
| 34 | knowledge of the CeCILL license and that you accept its terms.
|
|---|
| 35 |
|
|---|
| 36 | */
|
|---|
| 37 |
|
|---|
| 38 | #include <unistd.h>
|
|---|
| 39 | #include <stdlib.h>
|
|---|
| 40 | #include <string.h>
|
|---|
| 41 | #include <assert.h>
|
|---|
| 42 | #include <stdio.h>
|
|---|
| 43 | #include <time.h>
|
|---|
| 44 | #include <sys/types.h>
|
|---|
| 45 | #include <sys/time.h>
|
|---|
| 46 |
|
|---|
| 47 | #include "utils.h"
|
|---|
| 48 | // #include "parametres.h"
|
|---|
| 49 | real_t **
|
|---|
| 50 | allocate(int imin, int imax, int nvar) {
|
|---|
| 51 | int i;
|
|---|
| 52 |
|
|---|
| 53 | real_t **r = (real_t **) calloc(nvar, sizeof(real_t *));
|
|---|
| 54 | assert(r != NULL);
|
|---|
| 55 | for (i = 0; i < nvar; i++) {
|
|---|
| 56 | r[i] = DMalloc(imax - imin + 1 + MallocGuard);
|
|---|
| 57 | }
|
|---|
| 58 | return r;
|
|---|
| 59 | }
|
|---|
| 60 |
|
|---|
| 61 | #ifndef __MIC__
|
|---|
| 62 | #define NUMA_ALLOC 0
|
|---|
| 63 | #endif
|
|---|
| 64 |
|
|---|
| 65 | #ifdef __MIC__
|
|---|
| 66 | #define MEMSET 1
|
|---|
| 67 | #else
|
|---|
| 68 | #define MEMSET 0
|
|---|
| 69 | #endif
|
|---|
| 70 | #if NUMA_ALLOC==1
|
|---|
| 71 | #include <numa.h>
|
|---|
| 72 | #endif
|
|---|
| 73 |
|
|---|
| 74 | void DFree(real_t ** adr, size_t n)
|
|---|
| 75 | {
|
|---|
| 76 | #if NUMA_ALLOC == 1
|
|---|
| 77 | numa_free(*adr, sizeof(real_t) * (n + MallocGuard));
|
|---|
| 78 | #else
|
|---|
| 79 | free(*adr);
|
|---|
| 80 | #endif
|
|---|
| 81 | *adr = NULL;
|
|---|
| 82 | }
|
|---|
| 83 | void IFree(int ** adr, size_t n)
|
|---|
| 84 | {
|
|---|
| 85 | #if NUMA_ALLOC == 1
|
|---|
| 86 | numa_free(*adr, sizeof(int) * (n + MallocGuard));
|
|---|
| 87 | #else
|
|---|
| 88 | free(*adr);
|
|---|
| 89 | #endif
|
|---|
| 90 | *adr = NULL;
|
|---|
| 91 | }
|
|---|
| 92 |
|
|---|
| 93 | real_t *
|
|---|
| 94 | DMalloc(size_t n) {
|
|---|
| 95 | size_t i;
|
|---|
| 96 | #if NUMA_ALLOC == 1
|
|---|
| 97 | real_t *r = (real_t *) numa_alloc_interleaved((n + MallocGuard) * sizeof(real_t));
|
|---|
| 98 | #else
|
|---|
| 99 | real_t *r = (real_t *) calloc((n + MallocGuard), sizeof(real_t));
|
|---|
| 100 | #endif
|
|---|
| 101 | assert(r != NULL);
|
|---|
| 102 |
|
|---|
| 103 | #if MEMSET == 1
|
|---|
| 104 | memset(r, 1, n * sizeof(real_t));
|
|---|
| 105 | #else
|
|---|
| 106 | #ifndef NOTOUCHPAGE
|
|---|
| 107 | #pragma omp parallel for private(i) shared(r)
|
|---|
| 108 | for (i = 0; i < n; i++)
|
|---|
| 109 | r[i] = 0.0L;
|
|---|
| 110 | #endif
|
|---|
| 111 | #endif
|
|---|
| 112 | return r;
|
|---|
| 113 | }
|
|---|
| 114 |
|
|---|
| 115 | int *
|
|---|
| 116 | IMalloc(size_t n) {
|
|---|
| 117 | size_t i;
|
|---|
| 118 | #if NUMA_ALLOC == 1
|
|---|
| 119 | int *r = (int *) numa_alloc((n + MallocGuard) * sizeof(int));
|
|---|
| 120 | #else
|
|---|
| 121 | int *r = (int *) calloc((n + MallocGuard), sizeof(int));
|
|---|
| 122 | #endif
|
|---|
| 123 | assert(r != NULL);
|
|---|
| 124 |
|
|---|
| 125 | #if MEMSET == 1
|
|---|
| 126 | memset(r, 1, n * sizeof(int));
|
|---|
| 127 | #else
|
|---|
| 128 | #pragma omp parallel for private(i) shared(r)
|
|---|
| 129 | for (i = 0; i < n; i++)
|
|---|
| 130 | r[i] = 0;
|
|---|
| 131 | #endif
|
|---|
| 132 | return r;
|
|---|
| 133 | }
|
|---|
| 134 |
|
|---|
| 135 |
|
|---|
| 136 | #include "parametres.h"
|
|---|
| 137 | #define VALPERLINE 16
|
|---|
| 138 | void
|
|---|
| 139 | printuoldf(FILE * fic, const hydroparam_t H, hydrovar_t * Hv) {
|
|---|
| 140 | int i, j, nvar;
|
|---|
| 141 | for (nvar = 0; nvar < H.nvar; nvar++) {
|
|---|
| 142 | fprintf(fic, "=uold %d >\n", nvar);
|
|---|
| 143 | for (j = 0; j < H.nyt; j++) {
|
|---|
| 144 | int nbr = 1;
|
|---|
| 145 | for (i = 0; i < H.nxt; i++) {
|
|---|
| 146 | fprintf(fic, "%12.4e ", Hv->uold[IHv(i, j, nvar)]);
|
|---|
| 147 | nbr++;
|
|---|
| 148 | if (nbr == VALPERLINE) {
|
|---|
| 149 | fprintf(fic, "\n");
|
|---|
| 150 | fflush(fic);
|
|---|
| 151 | nbr = 1;
|
|---|
| 152 | }
|
|---|
| 153 | }
|
|---|
| 154 | if (nbr != 1)
|
|---|
| 155 | fprintf(fic, "\n");
|
|---|
| 156 | // fprintf(fic, "%%\n");
|
|---|
| 157 | fflush(fic);
|
|---|
| 158 | }
|
|---|
| 159 | }
|
|---|
| 160 | }
|
|---|
| 161 |
|
|---|
| 162 | void
|
|---|
| 163 | printarray(FILE * fic, real_t *a, int n, const char *nom, const hydroparam_t H) {
|
|---|
| 164 | real_t (*ptr)[H.nxyt] = (real_t (*)[H.nxyt]) a;
|
|---|
| 165 | long i, j, nbr = 1;
|
|---|
| 166 | fprintf(fic, "=%s >\n", nom);
|
|---|
| 167 | for (j = 0; j < H.nxystep; j++) {
|
|---|
| 168 | nbr = 1;
|
|---|
| 169 | for (i = 0; i < n; i++) {
|
|---|
| 170 | fprintf(fic, "%12.4e ", ptr[j][i]);
|
|---|
| 171 | nbr++;
|
|---|
| 172 | if (nbr == VALPERLINE) {
|
|---|
| 173 | fprintf(fic, "\n");
|
|---|
| 174 | nbr = 1;
|
|---|
| 175 | }
|
|---|
| 176 | }
|
|---|
| 177 | if (nbr != 1)
|
|---|
| 178 | fprintf(fic, "\n");
|
|---|
| 179 | }
|
|---|
| 180 | fprintf(fic, "\n");
|
|---|
| 181 | }
|
|---|
| 182 |
|
|---|
| 183 | void
|
|---|
| 184 | printarrayi(FILE * fic, int *a, int n, const char *nom) {
|
|---|
| 185 | int i, nbr = 1;
|
|---|
| 186 | fprintf(fic, "=%s >\n", nom);
|
|---|
| 187 | for (i = 0; i < n; i++) {
|
|---|
| 188 | fprintf(fic, "%4d ", a[i]);
|
|---|
| 189 | nbr++;
|
|---|
| 190 | if (nbr == VALPERLINE) {
|
|---|
| 191 | fprintf(fic, "\n");
|
|---|
| 192 | nbr = 1;
|
|---|
| 193 | }
|
|---|
| 194 | }
|
|---|
| 195 | if (nbr != 1)
|
|---|
| 196 | fprintf(fic, "\n");
|
|---|
| 197 | }
|
|---|
| 198 |
|
|---|
| 199 | void
|
|---|
| 200 | printarrayv(FILE * fic, real_t *a, int n, const char *nom, const hydroparam_t H) {
|
|---|
| 201 | int i, nbr = 1;
|
|---|
| 202 | int nvar;
|
|---|
| 203 | fprintf(fic, "=%s >\n", nom);
|
|---|
| 204 | real_t (*ptr)[H.nxyt] = (real_t (*)[H.nxyt]) a;
|
|---|
| 205 | for (nvar = 0; nvar < H.nvar; nvar++) {
|
|---|
| 206 | nbr = 1;
|
|---|
| 207 | for (i = 0; i < n; i++) {
|
|---|
| 208 | fprintf(fic, "%12.4e ", ptr[nvar][i]);
|
|---|
| 209 | nbr++;
|
|---|
| 210 | if (nbr == VALPERLINE) {
|
|---|
| 211 | fprintf(fic, "\n");
|
|---|
| 212 | nbr = 1;
|
|---|
| 213 | }
|
|---|
| 214 | }
|
|---|
| 215 | if (nbr != 1)
|
|---|
| 216 | fprintf(fic, "\n");
|
|---|
| 217 | fprintf(fic, "---\n");
|
|---|
| 218 | }
|
|---|
| 219 | }
|
|---|
| 220 |
|
|---|
| 221 | void
|
|---|
| 222 | printarrayv2(FILE * fic, real_t *a, int n, const char *nom, const hydroparam_t H) {
|
|---|
| 223 | int i, j, nbr = 1;
|
|---|
| 224 | int nvar;
|
|---|
| 225 | fprintf(fic, "=%s >\n#", nom);
|
|---|
| 226 | real_t (*ptr)[H.nxystep][H.nxyt] = (real_t (*)[H.nxystep][H.nxyt]) a;
|
|---|
| 227 | for (nvar = 0; nvar < H.nvar; nvar++) {
|
|---|
| 228 | for (j = 0; j < H.nxystep; j++) {
|
|---|
| 229 | nbr = 1;
|
|---|
| 230 | for (i = 0; i < n; i++) {
|
|---|
| 231 | fprintf(fic, "%12.4le ", ptr[nvar][j][i]);
|
|---|
| 232 | nbr++;
|
|---|
| 233 | if (nbr == VALPERLINE) {
|
|---|
| 234 | fprintf(fic, "\n#");
|
|---|
| 235 | nbr = 1;
|
|---|
| 236 | }
|
|---|
| 237 | }
|
|---|
| 238 | if (nbr != 1)
|
|---|
| 239 | fprintf(fic, "@\n#");
|
|---|
| 240 | }
|
|---|
| 241 | fprintf(fic, "-J-\n#");
|
|---|
| 242 | }
|
|---|
| 243 | fprintf(fic, "---\n");
|
|---|
| 244 | }
|
|---|
| 245 |
|
|---|
| 246 | void
|
|---|
| 247 | timeToString(char *buf, const double timeInS) {
|
|---|
| 248 | char ctenth[10];
|
|---|
| 249 | int hour = (int) (timeInS / 3600.0);
|
|---|
| 250 | int minute = (int) ((timeInS - hour * 3600) / 60.0);
|
|---|
| 251 | int second = (int) (timeInS - hour * 3600 - minute * 60);
|
|---|
| 252 | float tenth = (float) (timeInS - hour * 3600 - minute * 60 - second);
|
|---|
| 253 | sprintf(ctenth, "%.3f", tenth);
|
|---|
| 254 | sprintf(buf, "%02d:%02d:%02d%s", hour, minute, second, &ctenth[1]);
|
|---|
| 255 | }
|
|---|
| 256 |
|
|---|
| 257 |
|
|---|
| 258 | // double
|
|---|
| 259 | // cclock(void) {
|
|---|
| 260 | // const double micro = 1.0e-06; /* Conversion constant */
|
|---|
| 261 | // static long start = 0L, startu;
|
|---|
| 262 | // struct timeval tp; /* Structure used by gettimeofday */
|
|---|
| 263 | // double wall_time; /* To hold the result */
|
|---|
| 264 | // if (gettimeofday(&tp, NULL) == -1)
|
|---|
| 265 | // wall_time = -1.0e0;
|
|---|
| 266 |
|
|---|
| 267 | // else if (!start) {
|
|---|
| 268 | // start = tp.tv_sec;
|
|---|
| 269 | // startu = tp.tv_usec;
|
|---|
| 270 | // wall_time = 0.0e0;
|
|---|
| 271 | // } else
|
|---|
| 272 | // wall_time = (double) (tp.tv_sec - start) + micro * (tp.tv_usec - startu);
|
|---|
| 273 | // return wall_time;
|
|---|
| 274 | // }
|
|---|
| 275 |
|
|---|
| 276 |
|
|---|
| 277 | //EOF
|
|---|