1 2 /* 3 Code for manipulating distributed regular arrays in parallel. 4 */ 5 6 #include <petsc-private/daimpl.h> /*I "petscdmda.h" I*/ 7 8 /* 9 This allows the DMDA vectors to properly tell MATLAB their dimensions 10 */ 11 #if defined(PETSC_HAVE_MATLAB_ENGINE) 12 #include <engine.h> /* MATLAB include file */ 13 #include <mex.h> /* MATLAB include file */ 14 EXTERN_C_BEGIN 15 #undef __FUNCT__ 16 #define __FUNCT__ "VecMatlabEnginePut_DA2d" 17 PetscErrorCode VecMatlabEnginePut_DA2d(PetscObject obj,void *mengine) 18 { 19 PetscErrorCode ierr; 20 PetscInt n,m; 21 Vec vec = (Vec)obj; 22 PetscScalar *array; 23 mxArray *mat; 24 DM da; 25 26 PetscFunctionBegin; 27 ierr = PetscObjectQuery((PetscObject)vec,"DM",(PetscObject*)&da);CHKERRQ(ierr); 28 if (!da) SETERRQ(((PetscObject)vec)->comm,PETSC_ERR_ARG_WRONGSTATE,"Vector not associated with a DMDA"); 29 ierr = DMDAGetGhostCorners(da,0,0,0,&m,&n,0);CHKERRQ(ierr); 30 31 ierr = VecGetArray(vec,&array);CHKERRQ(ierr); 32 #if !defined(PETSC_USE_COMPLEX) 33 mat = mxCreateDoubleMatrix(m,n,mxREAL); 34 #else 35 mat = mxCreateDoubleMatrix(m,n,mxCOMPLEX); 36 #endif 37 ierr = PetscMemcpy(mxGetPr(mat),array,n*m*sizeof(PetscScalar));CHKERRQ(ierr); 38 ierr = PetscObjectName(obj);CHKERRQ(ierr); 39 engPutVariable((Engine *)mengine,obj->name,mat); 40 41 ierr = VecRestoreArray(vec,&array);CHKERRQ(ierr); 42 PetscFunctionReturn(0); 43 } 44 EXTERN_C_END 45 #endif 46 47 48 #undef __FUNCT__ 49 #define __FUNCT__ "DMCreateLocalVector_DA" 50 PetscErrorCode DMCreateLocalVector_DA(DM da,Vec* g) 51 { 52 PetscErrorCode ierr; 53 DM_DA *dd = (DM_DA*)da->data; 54 55 PetscFunctionBegin; 56 PetscValidHeaderSpecific(da,DM_CLASSID,1); 57 PetscValidPointer(g,2); 58 ierr = VecCreate(PETSC_COMM_SELF,g);CHKERRQ(ierr); 59 if (dd->defaultSection) { 60 PetscInt localSize; 61 62 ierr = PetscSectionGetStorageSize(dd->defaultSection, &localSize);CHKERRQ(ierr); 63 ierr = VecSetSizes(*g, localSize, PETSC_DETERMINE);CHKERRQ(ierr); 64 } else { 65 ierr = VecSetSizes(*g,dd->nlocal,PETSC_DETERMINE);CHKERRQ(ierr); 66 ierr = VecSetBlockSize(*g,dd->w);CHKERRQ(ierr); 67 } 68 ierr = VecSetType(*g,da->vectype);CHKERRQ(ierr); 69 ierr = PetscObjectCompose((PetscObject)*g,"DM",(PetscObject)da);CHKERRQ(ierr); 70 #if defined(PETSC_HAVE_MATLAB_ENGINE) 71 if (dd->w == 1 && dd->dim == 2) { 72 ierr = PetscObjectComposeFunctionDynamic((PetscObject)*g,"PetscMatlabEnginePut_C","VecMatlabEnginePut_DA2d",VecMatlabEnginePut_DA2d);CHKERRQ(ierr); 73 } 74 #endif 75 PetscFunctionReturn(0); 76 } 77 78 #undef __FUNCT__ 79 #define __FUNCT__ "DMDACreateSection" 80 /*@C 81 DMDACreateSection - Create a PetscSection inside the DMDA that describes data layout. This allows multiple fields with 82 different numbers of dofs on vertices, cells, and faces in each direction. 83 84 Input Parameters: 85 + dm- The DMDA 86 . numFields - The number of fields 87 . numComp - The number of components in each field, or PETSC_NULL for 1 88 . numVertexDof - The number of dofs per vertex for each field, or PETSC_NULL 89 . numFaceDof - The number of dofs per face for each field and direction, or PETSC_NULL 90 - numCellDof - The number of dofs per cell for each field, or PETSC_NULL 91 92 Level: developer 93 94 Note: 95 The default DMDA numbering is as follows: 96 97 - Cells: [0, nC) 98 - Vertices: [nC, nC+nV) 99 - X-Faces: [nC+nV, nC+nV+nXF) 100 - Y-Faces: [nC+nV+nXF, nC+nV+nXF+nYF) 101 - Z-Faces: [nC+nV+nXF+nYF, nC+nV+nXF+nYF+nZF) 102 103 We interpret the default DMDA partition as a cell partition, and the data assignment as a cell assignment. 104 @*/ 105 PetscErrorCode DMDACreateSection(DM dm, PetscInt numComp[], PetscInt numVertexDof[], PetscInt numFaceDof[], PetscInt numCellDof[]) 106 { 107 DM_DA *da = (DM_DA *) dm->data; 108 PetscInt mx = (da->Xe - da->Xs)/da->w, my = da->Ye - da->Ys, mz = da->Ze - da->Zs; 109 PetscInt nC = (mx )*(da->dim > 1 ? (my )*(da->dim > 2 ? (mz ) : 1) : 1); 110 PetscInt nV = (mx+1)*(da->dim > 1 ? (my+1)*(da->dim > 2 ? (mz+1) : 1) : 1); 111 PetscInt nXF = (mx+1)*(da->dim > 1 ? (my )*(da->dim > 2 ? (mz ) : 1) : 1); 112 PetscInt nYF = (mx )*(da->dim > 1 ? (my+1)*(da->dim > 2 ? (mz ) : 1) : 0); 113 PetscInt nZF = (mx )*(da->dim > 1 ? (my )*(da->dim > 2 ? (mz+1) : 0) : 0); 114 PetscInt cStart = 0, cEnd = cStart+nC; 115 PetscInt vStart = cEnd, vEnd = vStart+nV; 116 PetscInt xfStart = vEnd, xfEnd = xfStart+nXF; 117 PetscInt yfStart = xfEnd, yfEnd = yfStart+nYF; 118 PetscInt zfStart = yfEnd, zfEnd = zfStart+nZF; 119 PetscInt pStart = 0, pEnd = zfEnd; 120 PetscInt numFields, numVertexTotDof = 0, numCellTotDof = 0, numFaceTotDof[3] = {0, 0, 0}; 121 PetscInt f, v, c, xf, yf, zf; 122 PetscErrorCode ierr; 123 124 PetscFunctionBegin; 125 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 126 ierr = DMDAGetInfo(dm, 0,0,0,0,0,0,0, &numFields, 0,0,0,0,0);CHKERRQ(ierr); 127 for(f = 0; f < numFields; ++f) { 128 if (numVertexDof) {numVertexTotDof += numVertexDof[f];} 129 if (numCellDof) {numCellTotDof += numCellDof[f];} 130 if (numFaceDof) {numFaceTotDof[0] += numFaceDof[f*da->dim+0]; 131 numFaceTotDof[1] += da->dim > 1 ? numFaceDof[f*da->dim+1] : 0; 132 numFaceTotDof[2] += da->dim > 2 ? numFaceDof[f*da->dim+2] : 0;} 133 } 134 ierr = PetscSectionCreate(((PetscObject) dm)->comm, &da->defaultSection);CHKERRQ(ierr); 135 if (numFields > 1) { 136 ierr = PetscSectionSetNumFields(da->defaultSection, numFields);CHKERRQ(ierr); 137 for(f = 0; f < numFields; ++f) { 138 const char *name; 139 140 ierr = DMDAGetFieldName(dm, f, &name);CHKERRQ(ierr); 141 ierr = PetscSectionSetFieldName(da->defaultSection, f, name);CHKERRQ(ierr); 142 if (numComp) { 143 ierr = PetscSectionSetFieldComponents(da->defaultSection, f, numComp[f]);CHKERRQ(ierr); 144 } 145 } 146 } else { 147 numFields = 0; 148 } 149 ierr = PetscSectionSetChart(da->defaultSection, pStart, pEnd);CHKERRQ(ierr); 150 if (numVertexDof) { 151 for(v = vStart; v < vEnd; ++v) { 152 for(f = 0; f < numFields; ++f) { 153 ierr = PetscSectionSetFieldDof(da->defaultSection, v, f, numVertexDof[f]);CHKERRQ(ierr); 154 } 155 ierr = PetscSectionSetDof(da->defaultSection, v, numVertexTotDof);CHKERRQ(ierr); 156 } 157 } 158 if (numFaceDof) { 159 for(xf = xfStart; xf < xfEnd; ++xf) { 160 for(f = 0; f < numFields; ++f) { 161 ierr = PetscSectionSetFieldDof(da->defaultSection, xf, f, numFaceDof[f*da->dim+0]);CHKERRQ(ierr); 162 } 163 ierr = PetscSectionSetDof(da->defaultSection, xf, numFaceTotDof[0]);CHKERRQ(ierr); 164 } 165 for(yf = yfStart; yf < yfEnd; ++yf) { 166 for(f = 0; f < numFields; ++f) { 167 ierr = PetscSectionSetFieldDof(da->defaultSection, yf, f, numFaceDof[f*da->dim+1]);CHKERRQ(ierr); 168 } 169 ierr = PetscSectionSetDof(da->defaultSection, yf, numFaceTotDof[1]);CHKERRQ(ierr); 170 } 171 for(zf = zfStart; zf < zfEnd; ++zf) { 172 for(f = 0; f < numFields; ++f) { 173 ierr = PetscSectionSetFieldDof(da->defaultSection, zf, f, numFaceDof[f*da->dim+2]);CHKERRQ(ierr); 174 } 175 ierr = PetscSectionSetDof(da->defaultSection, zf, numFaceTotDof[2]);CHKERRQ(ierr); 176 } 177 } 178 if (numCellDof) { 179 for(c = cStart; c < cEnd; ++c) { 180 for(f = 0; f < numFields; ++f) { 181 ierr = PetscSectionSetFieldDof(da->defaultSection, c, f, numCellDof[f]);CHKERRQ(ierr); 182 } 183 ierr = PetscSectionSetDof(da->defaultSection, c, numCellTotDof);CHKERRQ(ierr); 184 } 185 } 186 ierr = PetscSectionSetUp(da->defaultSection);CHKERRQ(ierr); 187 da->defaultGlobalSection = da->defaultSection; 188 PetscFunctionReturn(0); 189 } 190 191 /* ------------------------------------------------------------------- */ 192 #if defined(PETSC_HAVE_ADIC) 193 194 EXTERN_C_BEGIN 195 #include <adic/ad_utils.h> 196 EXTERN_C_END 197 198 #undef __FUNCT__ 199 #define __FUNCT__ "DMDAGetAdicArray" 200 /*@C 201 DMDAGetAdicArray - Gets an array of derivative types for a DMDA 202 203 Input Parameter: 204 + da - information about my local patch 205 - ghosted - do you want arrays for the ghosted or nonghosted patch 206 207 Output Parameters: 208 + vptr - array data structured to be passed to ad_FormFunctionLocal() 209 . array_start - actual start of 1d array of all values that adiC can access directly (may be null) 210 - tdof - total number of degrees of freedom represented in array_start (may be null) 211 212 Notes: 213 The vector values are NOT initialized and may have garbage in them, so you may need 214 to zero them. 215 216 Returns the same type of object as the DMDAVecGetArray() except its elements are 217 derivative types instead of PetscScalars 218 219 Level: advanced 220 221 .seealso: DMDARestoreAdicArray() 222 223 @*/ 224 PetscErrorCode DMDAGetAdicArray(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 225 { 226 PetscErrorCode ierr; 227 PetscInt j,i,deriv_type_size,xs,ys,xm,ym,zs,zm,itdof; 228 char *iarray_start; 229 void **iptr = (void**)vptr; 230 DM_DA *dd = (DM_DA*)da->data; 231 232 PetscFunctionBegin; 233 PetscValidHeaderSpecific(da,DM_CLASSID,1); 234 if (ghosted) { 235 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 236 if (dd->adarrayghostedin[i]) { 237 *iptr = dd->adarrayghostedin[i]; 238 iarray_start = (char*)dd->adstartghostedin[i]; 239 itdof = dd->ghostedtdof; 240 dd->adarrayghostedin[i] = PETSC_NULL; 241 dd->adstartghostedin[i] = PETSC_NULL; 242 243 goto done; 244 } 245 } 246 xs = dd->Xs; 247 ys = dd->Ys; 248 zs = dd->Zs; 249 xm = dd->Xe-dd->Xs; 250 ym = dd->Ye-dd->Ys; 251 zm = dd->Ze-dd->Zs; 252 } else { 253 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 254 if (dd->adarrayin[i]) { 255 *iptr = dd->adarrayin[i]; 256 iarray_start = (char*)dd->adstartin[i]; 257 itdof = dd->tdof; 258 dd->adarrayin[i] = PETSC_NULL; 259 dd->adstartin[i] = PETSC_NULL; 260 261 goto done; 262 } 263 } 264 xs = dd->xs; 265 ys = dd->ys; 266 zs = dd->zs; 267 xm = dd->xe-dd->xs; 268 ym = dd->ye-dd->ys; 269 zm = dd->ze-dd->zs; 270 } 271 deriv_type_size = PetscADGetDerivTypeSize(); 272 273 switch (dd->dim) { 274 case 1: { 275 void *ptr; 276 itdof = xm; 277 278 ierr = PetscMalloc(xm*deriv_type_size,&iarray_start);CHKERRQ(ierr); 279 280 ptr = (void*)(iarray_start - xs*deriv_type_size); 281 *iptr = (void*)ptr; 282 break;} 283 case 2: { 284 void **ptr; 285 itdof = xm*ym; 286 287 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*deriv_type_size,&iarray_start);CHKERRQ(ierr); 288 289 ptr = (void**)(iarray_start + xm*ym*deriv_type_size - ys*sizeof(void*)); 290 for(j=ys;j<ys+ym;j++) { 291 ptr[j] = iarray_start + deriv_type_size*(xm*(j-ys) - xs); 292 } 293 *iptr = (void*)ptr; 294 break;} 295 case 3: { 296 void ***ptr,**bptr; 297 itdof = xm*ym*zm; 298 299 ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*deriv_type_size,&iarray_start);CHKERRQ(ierr); 300 301 ptr = (void***)(iarray_start + xm*ym*zm*deriv_type_size - zs*sizeof(void*)); 302 bptr = (void**)(iarray_start + xm*ym*zm*deriv_type_size + zm*sizeof(void**)); 303 304 for(i=zs;i<zs+zm;i++) { 305 ptr[i] = bptr + ((i-zs)*ym - ys); 306 } 307 for (i=zs; i<zs+zm; i++) { 308 for (j=ys; j<ys+ym; j++) { 309 ptr[i][j] = iarray_start + deriv_type_size*(xm*ym*(i-zs) + xm*(j-ys) - xs); 310 } 311 } 312 313 *iptr = (void*)ptr; 314 break;} 315 default: 316 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 317 } 318 319 done: 320 /* add arrays to the checked out list */ 321 if (ghosted) { 322 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 323 if (!dd->adarrayghostedout[i]) { 324 dd->adarrayghostedout[i] = *iptr ; 325 dd->adstartghostedout[i] = iarray_start; 326 dd->ghostedtdof = itdof; 327 break; 328 } 329 } 330 } else { 331 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 332 if (!dd->adarrayout[i]) { 333 dd->adarrayout[i] = *iptr ; 334 dd->adstartout[i] = iarray_start; 335 dd->tdof = itdof; 336 break; 337 } 338 } 339 } 340 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Too many DMDA ADIC arrays obtained"); 341 if (tdof) *tdof = itdof; 342 if (array_start) *(void**)array_start = iarray_start; 343 PetscFunctionReturn(0); 344 } 345 346 #undef __FUNCT__ 347 #define __FUNCT__ "DMDARestoreAdicArray" 348 /*@C 349 DMDARestoreAdicArray - Restores an array of derivative types for a DMDA 350 351 Input Parameter: 352 + da - information about my local patch 353 - ghosted - do you want arrays for the ghosted or nonghosted patch 354 355 Output Parameters: 356 + ptr - array data structured to be passed to ad_FormFunctionLocal() 357 . array_start - actual start of 1d array of all values that adiC can access directly 358 - tdof - total number of degrees of freedom represented in array_start 359 360 Level: advanced 361 362 .seealso: DMDAGetAdicArray() 363 364 @*/ 365 PetscErrorCode DMDARestoreAdicArray(DM da,PetscBool ghosted,void *ptr,void *array_start,PetscInt *tdof) 366 { 367 PetscInt i; 368 void **iptr = (void**)ptr,iarray_start = 0; 369 370 PetscFunctionBegin; 371 PetscValidHeaderSpecific(da,DM_CLASSID,1); 372 if (ghosted) { 373 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 374 if (dd->adarrayghostedout[i] == *iptr) { 375 iarray_start = dd->adstartghostedout[i]; 376 dd->adarrayghostedout[i] = PETSC_NULL; 377 dd->adstartghostedout[i] = PETSC_NULL; 378 break; 379 } 380 } 381 if (!iarray_start) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not find array in checkout list"); 382 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 383 if (!dd->adarrayghostedin[i]){ 384 dd->adarrayghostedin[i] = *iptr; 385 dd->adstartghostedin[i] = iarray_start; 386 break; 387 } 388 } 389 } else { 390 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 391 if (dd->adarrayout[i] == *iptr) { 392 iarray_start = dd->adstartout[i]; 393 dd->adarrayout[i] = PETSC_NULL; 394 dd->adstartout[i] = PETSC_NULL; 395 break; 396 } 397 } 398 if (!iarray_start) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not find array in checkout list"); 399 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 400 if (!dd->adarrayin[i]){ 401 dd->adarrayin[i] = *iptr; 402 dd->adstartin[i] = iarray_start; 403 break; 404 } 405 } 406 } 407 PetscFunctionReturn(0); 408 } 409 410 #undef __FUNCT__ 411 #define __FUNCT__ "ad_DAGetArray" 412 PetscErrorCode ad_DAGetArray(DM da,PetscBool ghosted,void *iptr) 413 { 414 PetscErrorCode ierr; 415 PetscFunctionBegin; 416 ierr = DMDAGetAdicArray(da,ghosted,iptr,0,0);CHKERRQ(ierr); 417 PetscFunctionReturn(0); 418 } 419 420 #undef __FUNCT__ 421 #define __FUNCT__ "ad_DARestoreArray" 422 PetscErrorCode ad_DARestoreArray(DM da,PetscBool ghosted,void *iptr) 423 { 424 PetscErrorCode ierr; 425 PetscFunctionBegin; 426 ierr = DMDARestoreAdicArray(da,ghosted,iptr,0,0);CHKERRQ(ierr); 427 PetscFunctionReturn(0); 428 } 429 430 #endif 431 432 #undef __FUNCT__ 433 #define __FUNCT__ "DMDAGetArray" 434 /*@C 435 DMDAGetArray - Gets a work array for a DMDA 436 437 Input Parameter: 438 + da - information about my local patch 439 - ghosted - do you want arrays for the ghosted or nonghosted patch 440 441 Output Parameters: 442 . vptr - array data structured 443 444 Note: The vector values are NOT initialized and may have garbage in them, so you may need 445 to zero them. 446 447 Level: advanced 448 449 .seealso: DMDARestoreArray(), DMDAGetAdicArray() 450 451 @*/ 452 PetscErrorCode DMDAGetArray(DM da,PetscBool ghosted,void *vptr) 453 { 454 PetscErrorCode ierr; 455 PetscInt j,i,xs,ys,xm,ym,zs,zm; 456 char *iarray_start; 457 void **iptr = (void**)vptr; 458 DM_DA *dd = (DM_DA*)da->data; 459 460 PetscFunctionBegin; 461 PetscValidHeaderSpecific(da,DM_CLASSID,1); 462 if (ghosted) { 463 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 464 if (dd->arrayghostedin[i]) { 465 *iptr = dd->arrayghostedin[i]; 466 iarray_start = (char*)dd->startghostedin[i]; 467 dd->arrayghostedin[i] = PETSC_NULL; 468 dd->startghostedin[i] = PETSC_NULL; 469 470 goto done; 471 } 472 } 473 xs = dd->Xs; 474 ys = dd->Ys; 475 zs = dd->Zs; 476 xm = dd->Xe-dd->Xs; 477 ym = dd->Ye-dd->Ys; 478 zm = dd->Ze-dd->Zs; 479 } else { 480 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 481 if (dd->arrayin[i]) { 482 *iptr = dd->arrayin[i]; 483 iarray_start = (char*)dd->startin[i]; 484 dd->arrayin[i] = PETSC_NULL; 485 dd->startin[i] = PETSC_NULL; 486 487 goto done; 488 } 489 } 490 xs = dd->xs; 491 ys = dd->ys; 492 zs = dd->zs; 493 xm = dd->xe-dd->xs; 494 ym = dd->ye-dd->ys; 495 zm = dd->ze-dd->zs; 496 } 497 498 switch (dd->dim) { 499 case 1: { 500 void *ptr; 501 502 ierr = PetscMalloc(xm*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 503 504 ptr = (void*)(iarray_start - xs*sizeof(PetscScalar)); 505 *iptr = (void*)ptr; 506 break;} 507 case 2: { 508 void **ptr; 509 510 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 511 512 ptr = (void**)(iarray_start + xm*ym*sizeof(PetscScalar) - ys*sizeof(void*)); 513 for(j=ys;j<ys+ym;j++) { 514 ptr[j] = iarray_start + sizeof(PetscScalar)*(xm*(j-ys) - xs); 515 } 516 *iptr = (void*)ptr; 517 break;} 518 case 3: { 519 void ***ptr,**bptr; 520 521 ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 522 523 ptr = (void***)(iarray_start + xm*ym*zm*sizeof(PetscScalar) - zs*sizeof(void*)); 524 bptr = (void**)(iarray_start + xm*ym*zm*sizeof(PetscScalar) + zm*sizeof(void**)); 525 for(i=zs;i<zs+zm;i++) { 526 ptr[i] = bptr + ((i-zs)*ym - ys); 527 } 528 for (i=zs; i<zs+zm; i++) { 529 for (j=ys; j<ys+ym; j++) { 530 ptr[i][j] = iarray_start + sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs); 531 } 532 } 533 534 *iptr = (void*)ptr; 535 break;} 536 default: 537 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 538 } 539 540 done: 541 /* add arrays to the checked out list */ 542 if (ghosted) { 543 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 544 if (!dd->arrayghostedout[i]) { 545 dd->arrayghostedout[i] = *iptr ; 546 dd->startghostedout[i] = iarray_start; 547 break; 548 } 549 } 550 } else { 551 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 552 if (!dd->arrayout[i]) { 553 dd->arrayout[i] = *iptr ; 554 dd->startout[i] = iarray_start; 555 break; 556 } 557 } 558 } 559 PetscFunctionReturn(0); 560 } 561 562 #undef __FUNCT__ 563 #define __FUNCT__ "DMDARestoreArray" 564 /*@C 565 DMDARestoreArray - Restores an array of derivative types for a DMDA 566 567 Input Parameter: 568 + da - information about my local patch 569 . ghosted - do you want arrays for the ghosted or nonghosted patch 570 - vptr - array data structured to be passed to ad_FormFunctionLocal() 571 572 Level: advanced 573 574 .seealso: DMDAGetArray(), DMDAGetAdicArray() 575 576 @*/ 577 PetscErrorCode DMDARestoreArray(DM da,PetscBool ghosted,void *vptr) 578 { 579 PetscInt i; 580 void **iptr = (void**)vptr,*iarray_start = 0; 581 DM_DA *dd = (DM_DA*)da->data; 582 583 PetscFunctionBegin; 584 PetscValidHeaderSpecific(da,DM_CLASSID,1); 585 if (ghosted) { 586 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 587 if (dd->arrayghostedout[i] == *iptr) { 588 iarray_start = dd->startghostedout[i]; 589 dd->arrayghostedout[i] = PETSC_NULL; 590 dd->startghostedout[i] = PETSC_NULL; 591 break; 592 } 593 } 594 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 595 if (!dd->arrayghostedin[i]){ 596 dd->arrayghostedin[i] = *iptr; 597 dd->startghostedin[i] = iarray_start; 598 break; 599 } 600 } 601 } else { 602 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 603 if (dd->arrayout[i] == *iptr) { 604 iarray_start = dd->startout[i]; 605 dd->arrayout[i] = PETSC_NULL; 606 dd->startout[i] = PETSC_NULL; 607 break; 608 } 609 } 610 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 611 if (!dd->arrayin[i]){ 612 dd->arrayin[i] = *iptr; 613 dd->startin[i] = iarray_start; 614 break; 615 } 616 } 617 } 618 PetscFunctionReturn(0); 619 } 620 621 #undef __FUNCT__ 622 #define __FUNCT__ "DMDAGetAdicMFArray" 623 /*@C 624 DMDAGetAdicMFArray - Gets an array of derivative types for a DMDA for matrix-free ADIC. 625 626 Input Parameter: 627 + da - information about my local patch 628 - ghosted - do you want arrays for the ghosted or nonghosted patch? 629 630 Output Parameters: 631 + vptr - array data structured to be passed to ad_FormFunctionLocal() 632 . array_start - actual start of 1d array of all values that adiC can access directly (may be null) 633 - tdof - total number of degrees of freedom represented in array_start (may be null) 634 635 Notes: 636 The vector values are NOT initialized and may have garbage in them, so you may need 637 to zero them. 638 639 This routine returns the same type of object as the DMDAVecGetArray(), except its 640 elements are derivative types instead of PetscScalars. 641 642 Level: advanced 643 644 .seealso: DMDARestoreAdicMFArray(), DMDAGetArray(), DMDAGetAdicArray() 645 646 @*/ 647 PetscErrorCode DMDAGetAdicMFArray(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 648 { 649 PetscErrorCode ierr; 650 PetscInt j,i,xs,ys,xm,ym,zs,zm,itdof = 0; 651 char *iarray_start; 652 void **iptr = (void**)vptr; 653 DM_DA *dd = (DM_DA*)da->data; 654 655 PetscFunctionBegin; 656 PetscValidHeaderSpecific(da,DM_CLASSID,1); 657 if (ghosted) { 658 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 659 if (dd->admfarrayghostedin[i]) { 660 *iptr = dd->admfarrayghostedin[i]; 661 iarray_start = (char*)dd->admfstartghostedin[i]; 662 itdof = dd->ghostedtdof; 663 dd->admfarrayghostedin[i] = PETSC_NULL; 664 dd->admfstartghostedin[i] = PETSC_NULL; 665 666 goto done; 667 } 668 } 669 xs = dd->Xs; 670 ys = dd->Ys; 671 zs = dd->Zs; 672 xm = dd->Xe-dd->Xs; 673 ym = dd->Ye-dd->Ys; 674 zm = dd->Ze-dd->Zs; 675 } else { 676 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 677 if (dd->admfarrayin[i]) { 678 *iptr = dd->admfarrayin[i]; 679 iarray_start = (char*)dd->admfstartin[i]; 680 itdof = dd->tdof; 681 dd->admfarrayin[i] = PETSC_NULL; 682 dd->admfstartin[i] = PETSC_NULL; 683 684 goto done; 685 } 686 } 687 xs = dd->xs; 688 ys = dd->ys; 689 zs = dd->zs; 690 xm = dd->xe-dd->xs; 691 ym = dd->ye-dd->ys; 692 zm = dd->ze-dd->zs; 693 } 694 695 switch (dd->dim) { 696 case 1: { 697 void *ptr; 698 itdof = xm; 699 700 ierr = PetscMalloc(xm*2*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 701 702 ptr = (void*)(iarray_start - xs*2*sizeof(PetscScalar)); 703 *iptr = (void*)ptr; 704 break;} 705 case 2: { 706 void **ptr; 707 itdof = xm*ym; 708 709 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*2*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 710 711 ptr = (void**)(iarray_start + xm*ym*2*sizeof(PetscScalar) - ys*sizeof(void*)); 712 for(j=ys;j<ys+ym;j++) { 713 ptr[j] = iarray_start + 2*sizeof(PetscScalar)*(xm*(j-ys) - xs); 714 } 715 *iptr = (void*)ptr; 716 break;} 717 case 3: { 718 void ***ptr,**bptr; 719 itdof = xm*ym*zm; 720 721 ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*2*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 722 723 ptr = (void***)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) - zs*sizeof(void*)); 724 bptr = (void**)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) + zm*sizeof(void**)); 725 for(i=zs;i<zs+zm;i++) { 726 ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*); 727 } 728 for (i=zs; i<zs+zm; i++) { 729 for (j=ys; j<ys+ym; j++) { 730 ptr[i][j] = iarray_start + 2*sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs); 731 } 732 } 733 734 *iptr = (void*)ptr; 735 break;} 736 default: 737 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 738 } 739 740 done: 741 /* add arrays to the checked out list */ 742 if (ghosted) { 743 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 744 if (!dd->admfarrayghostedout[i]) { 745 dd->admfarrayghostedout[i] = *iptr ; 746 dd->admfstartghostedout[i] = iarray_start; 747 dd->ghostedtdof = itdof; 748 break; 749 } 750 } 751 } else { 752 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 753 if (!dd->admfarrayout[i]) { 754 dd->admfarrayout[i] = *iptr ; 755 dd->admfstartout[i] = iarray_start; 756 dd->tdof = itdof; 757 break; 758 } 759 } 760 } 761 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(((PetscObject)da)->comm,PETSC_ERR_ARG_WRONG,"Too many DMDA ADIC arrays obtained"); 762 if (tdof) *tdof = itdof; 763 if (array_start) *(void**)array_start = iarray_start; 764 PetscFunctionReturn(0); 765 } 766 767 #undef __FUNCT__ 768 #define __FUNCT__ "DMDAGetAdicMFArray4" 769 PetscErrorCode DMDAGetAdicMFArray4(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 770 { 771 PetscErrorCode ierr; 772 PetscInt j,i,xs,ys,xm,ym,itdof = 0; 773 char *iarray_start; 774 void **iptr = (void**)vptr; 775 DM_DA *dd = (DM_DA*)da->data; 776 777 PetscFunctionBegin; 778 PetscValidHeaderSpecific(da,DM_CLASSID,1); 779 if (ghosted) { 780 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 781 if (dd->admfarrayghostedin[i]) { 782 *iptr = dd->admfarrayghostedin[i]; 783 iarray_start = (char*)dd->admfstartghostedin[i]; 784 itdof = dd->ghostedtdof; 785 dd->admfarrayghostedin[i] = PETSC_NULL; 786 dd->admfstartghostedin[i] = PETSC_NULL; 787 788 goto done; 789 } 790 } 791 xs = dd->Xs; 792 ys = dd->Ys; 793 xm = dd->Xe-dd->Xs; 794 ym = dd->Ye-dd->Ys; 795 } else { 796 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 797 if (dd->admfarrayin[i]) { 798 *iptr = dd->admfarrayin[i]; 799 iarray_start = (char*)dd->admfstartin[i]; 800 itdof = dd->tdof; 801 dd->admfarrayin[i] = PETSC_NULL; 802 dd->admfstartin[i] = PETSC_NULL; 803 804 goto done; 805 } 806 } 807 xs = dd->xs; 808 ys = dd->ys; 809 xm = dd->xe-dd->xs; 810 ym = dd->ye-dd->ys; 811 } 812 813 switch (dd->dim) { 814 case 2: { 815 void **ptr; 816 itdof = xm*ym; 817 818 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*5*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 819 820 ptr = (void**)(iarray_start + xm*ym*5*sizeof(PetscScalar) - ys*sizeof(void*)); 821 for(j=ys;j<ys+ym;j++) { 822 ptr[j] = iarray_start + 5*sizeof(PetscScalar)*(xm*(j-ys) - xs); 823 } 824 *iptr = (void*)ptr; 825 break;} 826 default: 827 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 828 } 829 830 done: 831 /* add arrays to the checked out list */ 832 if (ghosted) { 833 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 834 if (!dd->admfarrayghostedout[i]) { 835 dd->admfarrayghostedout[i] = *iptr ; 836 dd->admfstartghostedout[i] = iarray_start; 837 dd->ghostedtdof = itdof; 838 break; 839 } 840 } 841 } else { 842 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 843 if (!dd->admfarrayout[i]) { 844 dd->admfarrayout[i] = *iptr ; 845 dd->admfstartout[i] = iarray_start; 846 dd->tdof = itdof; 847 break; 848 } 849 } 850 } 851 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Too many DMDA ADIC arrays obtained"); 852 if (tdof) *tdof = itdof; 853 if (array_start) *(void**)array_start = iarray_start; 854 PetscFunctionReturn(0); 855 } 856 857 #undef __FUNCT__ 858 #define __FUNCT__ "DMDAGetAdicMFArray9" 859 PetscErrorCode DMDAGetAdicMFArray9(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 860 { 861 PetscErrorCode ierr; 862 PetscInt j,i,xs,ys,xm,ym,itdof = 0; 863 char *iarray_start; 864 void **iptr = (void**)vptr; 865 DM_DA *dd = (DM_DA*)da->data; 866 867 PetscFunctionBegin; 868 PetscValidHeaderSpecific(da,DM_CLASSID,1); 869 if (ghosted) { 870 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 871 if (dd->admfarrayghostedin[i]) { 872 *iptr = dd->admfarrayghostedin[i]; 873 iarray_start = (char*)dd->admfstartghostedin[i]; 874 itdof = dd->ghostedtdof; 875 dd->admfarrayghostedin[i] = PETSC_NULL; 876 dd->admfstartghostedin[i] = PETSC_NULL; 877 878 goto done; 879 } 880 } 881 xs = dd->Xs; 882 ys = dd->Ys; 883 xm = dd->Xe-dd->Xs; 884 ym = dd->Ye-dd->Ys; 885 } else { 886 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 887 if (dd->admfarrayin[i]) { 888 *iptr = dd->admfarrayin[i]; 889 iarray_start = (char*)dd->admfstartin[i]; 890 itdof = dd->tdof; 891 dd->admfarrayin[i] = PETSC_NULL; 892 dd->admfstartin[i] = PETSC_NULL; 893 894 goto done; 895 } 896 } 897 xs = dd->xs; 898 ys = dd->ys; 899 xm = dd->xe-dd->xs; 900 ym = dd->ye-dd->ys; 901 } 902 903 switch (dd->dim) { 904 case 2: { 905 void **ptr; 906 itdof = xm*ym; 907 908 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*10*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 909 910 ptr = (void**)(iarray_start + xm*ym*10*sizeof(PetscScalar) - ys*sizeof(void*)); 911 for(j=ys;j<ys+ym;j++) { 912 ptr[j] = iarray_start + 10*sizeof(PetscScalar)*(xm*(j-ys) - xs); 913 } 914 *iptr = (void*)ptr; 915 break;} 916 default: 917 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 918 } 919 920 done: 921 /* add arrays to the checked out list */ 922 if (ghosted) { 923 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 924 if (!dd->admfarrayghostedout[i]) { 925 dd->admfarrayghostedout[i] = *iptr ; 926 dd->admfstartghostedout[i] = iarray_start; 927 dd->ghostedtdof = itdof; 928 break; 929 } 930 } 931 } else { 932 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 933 if (!dd->admfarrayout[i]) { 934 dd->admfarrayout[i] = *iptr ; 935 dd->admfstartout[i] = iarray_start; 936 dd->tdof = itdof; 937 break; 938 } 939 } 940 } 941 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Too many DMDA ADIC arrays obtained"); 942 if (tdof) *tdof = itdof; 943 if (array_start) *(void**)array_start = iarray_start; 944 PetscFunctionReturn(0); 945 } 946 947 #undef __FUNCT__ 948 #define __FUNCT__ "DMDAGetAdicMFArrayb" 949 /*@C 950 DMDAGetAdicMFArrayb - Gets an array of derivative types for a DMDA for matrix-free ADIC. 951 952 Input Parameter: 953 + da - information about my local patch 954 - ghosted - do you want arrays for the ghosted or nonghosted patch? 955 956 Output Parameters: 957 + vptr - array data structured to be passed to ad_FormFunctionLocal() 958 . array_start - actual start of 1d array of all values that adiC can access directly (may be null) 959 - tdof - total number of degrees of freedom represented in array_start (may be null) 960 961 Notes: 962 The vector values are NOT initialized and may have garbage in them, so you may need 963 to zero them. 964 965 This routine returns the same type of object as the DMDAVecGetArray(), except its 966 elements are derivative types instead of PetscScalars. 967 968 Level: advanced 969 970 .seealso: DMDARestoreAdicMFArray(), DMDAGetArray(), DMDAGetAdicArray() 971 972 @*/ 973 PetscErrorCode DMDAGetAdicMFArrayb(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 974 { 975 PetscErrorCode ierr; 976 PetscInt j,i,xs,ys,xm,ym,zs,zm,itdof = 0; 977 char *iarray_start; 978 void **iptr = (void**)vptr; 979 DM_DA *dd = (DM_DA*)da->data; 980 PetscInt bs = dd->w,bs1 = bs+1; 981 982 PetscFunctionBegin; 983 PetscValidHeaderSpecific(da,DM_CLASSID,1); 984 if (ghosted) { 985 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 986 if (dd->admfarrayghostedin[i]) { 987 *iptr = dd->admfarrayghostedin[i]; 988 iarray_start = (char*)dd->admfstartghostedin[i]; 989 itdof = dd->ghostedtdof; 990 dd->admfarrayghostedin[i] = PETSC_NULL; 991 dd->admfstartghostedin[i] = PETSC_NULL; 992 993 goto done; 994 } 995 } 996 xs = dd->Xs; 997 ys = dd->Ys; 998 zs = dd->Zs; 999 xm = dd->Xe-dd->Xs; 1000 ym = dd->Ye-dd->Ys; 1001 zm = dd->Ze-dd->Zs; 1002 } else { 1003 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1004 if (dd->admfarrayin[i]) { 1005 *iptr = dd->admfarrayin[i]; 1006 iarray_start = (char*)dd->admfstartin[i]; 1007 itdof = dd->tdof; 1008 dd->admfarrayin[i] = PETSC_NULL; 1009 dd->admfstartin[i] = PETSC_NULL; 1010 1011 goto done; 1012 } 1013 } 1014 xs = dd->xs; 1015 ys = dd->ys; 1016 zs = dd->zs; 1017 xm = dd->xe-dd->xs; 1018 ym = dd->ye-dd->ys; 1019 zm = dd->ze-dd->zs; 1020 } 1021 1022 switch (dd->dim) { 1023 case 1: { 1024 void *ptr; 1025 itdof = xm; 1026 1027 ierr = PetscMalloc(xm*bs1*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 1028 1029 ptr = (void*)(iarray_start - xs*bs1*sizeof(PetscScalar)); 1030 *iptr = (void*)ptr; 1031 break;} 1032 case 2: { 1033 void **ptr; 1034 itdof = xm*ym; 1035 1036 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*bs1*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 1037 1038 ptr = (void**)(iarray_start + xm*ym*bs1*sizeof(PetscScalar) - ys*sizeof(void*)); 1039 for(j=ys;j<ys+ym;j++) { 1040 ptr[j] = iarray_start + bs1*sizeof(PetscScalar)*(xm*(j-ys) - xs); 1041 } 1042 *iptr = (void*)ptr; 1043 break;} 1044 case 3: { 1045 void ***ptr,**bptr; 1046 itdof = xm*ym*zm; 1047 1048 ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*bs1*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 1049 1050 ptr = (void***)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) - zs*sizeof(void*)); 1051 bptr = (void**)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) + zm*sizeof(void**)); 1052 for(i=zs;i<zs+zm;i++) { 1053 ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*); 1054 } 1055 for (i=zs; i<zs+zm; i++) { 1056 for (j=ys; j<ys+ym; j++) { 1057 ptr[i][j] = iarray_start + bs1*sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs); 1058 } 1059 } 1060 1061 *iptr = (void*)ptr; 1062 break;} 1063 default: 1064 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 1065 } 1066 1067 done: 1068 /* add arrays to the checked out list */ 1069 if (ghosted) { 1070 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1071 if (!dd->admfarrayghostedout[i]) { 1072 dd->admfarrayghostedout[i] = *iptr ; 1073 dd->admfstartghostedout[i] = iarray_start; 1074 dd->ghostedtdof = itdof; 1075 break; 1076 } 1077 } 1078 } else { 1079 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1080 if (!dd->admfarrayout[i]) { 1081 dd->admfarrayout[i] = *iptr ; 1082 dd->admfstartout[i] = iarray_start; 1083 dd->tdof = itdof; 1084 break; 1085 } 1086 } 1087 } 1088 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Too many DMDA ADIC arrays obtained"); 1089 if (tdof) *tdof = itdof; 1090 if (array_start) *(void**)array_start = iarray_start; 1091 PetscFunctionReturn(0); 1092 } 1093 1094 #undef __FUNCT__ 1095 #define __FUNCT__ "DMDARestoreAdicMFArray" 1096 /*@C 1097 DMDARestoreAdicMFArray - Restores an array of derivative types for a DMDA. 1098 1099 Input Parameter: 1100 + da - information about my local patch 1101 - ghosted - do you want arrays for the ghosted or nonghosted patch? 1102 1103 Output Parameters: 1104 + ptr - array data structure to be passed to ad_FormFunctionLocal() 1105 . array_start - actual start of 1d array of all values that adiC can access directly 1106 - tdof - total number of degrees of freedom represented in array_start 1107 1108 Level: advanced 1109 1110 .seealso: DMDAGetAdicArray() 1111 1112 @*/ 1113 PetscErrorCode DMDARestoreAdicMFArray(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 1114 { 1115 PetscInt i; 1116 void **iptr = (void**)vptr,*iarray_start = 0; 1117 DM_DA *dd = (DM_DA*)da->data; 1118 1119 PetscFunctionBegin; 1120 PetscValidHeaderSpecific(da,DM_CLASSID,1); 1121 if (ghosted) { 1122 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1123 if (dd->admfarrayghostedout[i] == *iptr) { 1124 iarray_start = dd->admfstartghostedout[i]; 1125 dd->admfarrayghostedout[i] = PETSC_NULL; 1126 dd->admfstartghostedout[i] = PETSC_NULL; 1127 break; 1128 } 1129 } 1130 if (!iarray_start) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not find array in checkout list"); 1131 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1132 if (!dd->admfarrayghostedin[i]){ 1133 dd->admfarrayghostedin[i] = *iptr; 1134 dd->admfstartghostedin[i] = iarray_start; 1135 break; 1136 } 1137 } 1138 } else { 1139 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1140 if (dd->admfarrayout[i] == *iptr) { 1141 iarray_start = dd->admfstartout[i]; 1142 dd->admfarrayout[i] = PETSC_NULL; 1143 dd->admfstartout[i] = PETSC_NULL; 1144 break; 1145 } 1146 } 1147 if (!iarray_start) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not find array in checkout list"); 1148 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1149 if (!dd->admfarrayin[i]){ 1150 dd->admfarrayin[i] = *iptr; 1151 dd->admfstartin[i] = iarray_start; 1152 break; 1153 } 1154 } 1155 } 1156 PetscFunctionReturn(0); 1157 } 1158 1159 #undef __FUNCT__ 1160 #define __FUNCT__ "admf_DAGetArray" 1161 PetscErrorCode admf_DAGetArray(DM da,PetscBool ghosted,void *iptr) 1162 { 1163 PetscErrorCode ierr; 1164 PetscFunctionBegin; 1165 ierr = DMDAGetAdicMFArray(da,ghosted,iptr,0,0);CHKERRQ(ierr); 1166 PetscFunctionReturn(0); 1167 } 1168 1169 #undef __FUNCT__ 1170 #define __FUNCT__ "admf_DARestoreArray" 1171 PetscErrorCode admf_DARestoreArray(DM da,PetscBool ghosted,void *iptr) 1172 { 1173 PetscErrorCode ierr; 1174 PetscFunctionBegin; 1175 ierr = DMDARestoreAdicMFArray(da,ghosted,iptr,0,0);CHKERRQ(ierr); 1176 PetscFunctionReturn(0); 1177 } 1178 1179