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 ierr = VecSetSizes(*g,dd->nlocal,PETSC_DETERMINE);CHKERRQ(ierr); 60 ierr = VecSetType(*g,da->vectype);CHKERRQ(ierr); 61 ierr = VecSetBlockSize(*g,dd->w);CHKERRQ(ierr); 62 ierr = PetscObjectCompose((PetscObject)*g,"DM",(PetscObject)da);CHKERRQ(ierr); 63 #if defined(PETSC_HAVE_MATLAB_ENGINE) 64 if (dd->w == 1 && dd->dim == 2) { 65 ierr = PetscObjectComposeFunctionDynamic((PetscObject)*g,"PetscMatlabEnginePut_C","VecMatlabEnginePut_DA2d",VecMatlabEnginePut_DA2d);CHKERRQ(ierr); 66 } 67 #endif 68 PetscFunctionReturn(0); 69 } 70 71 /* ------------------------------------------------------------------- */ 72 #if defined(PETSC_HAVE_ADIC) 73 74 EXTERN_C_BEGIN 75 #include <adic/ad_utils.h> 76 EXTERN_C_END 77 78 #undef __FUNCT__ 79 #define __FUNCT__ "DMDAGetAdicArray" 80 /*@C 81 DMDAGetAdicArray - Gets an array of derivative types for a DMDA 82 83 Input Parameter: 84 + da - information about my local patch 85 - ghosted - do you want arrays for the ghosted or nonghosted patch 86 87 Output Parameters: 88 + vptr - array data structured to be passed to ad_FormFunctionLocal() 89 . array_start - actual start of 1d array of all values that adiC can access directly (may be null) 90 - tdof - total number of degrees of freedom represented in array_start (may be null) 91 92 Notes: 93 The vector values are NOT initialized and may have garbage in them, so you may need 94 to zero them. 95 96 Returns the same type of object as the DMDAVecGetArray() except its elements are 97 derivative types instead of PetscScalars 98 99 Level: advanced 100 101 .seealso: DMDARestoreAdicArray() 102 103 @*/ 104 PetscErrorCode DMDAGetAdicArray(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 105 { 106 PetscErrorCode ierr; 107 PetscInt j,i,deriv_type_size,xs,ys,xm,ym,zs,zm,itdof; 108 char *iarray_start; 109 void **iptr = (void**)vptr; 110 DM_DA *dd = (DM_DA*)da->data; 111 112 PetscFunctionBegin; 113 PetscValidHeaderSpecific(da,DM_CLASSID,1); 114 if (ghosted) { 115 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 116 if (dd->adarrayghostedin[i]) { 117 *iptr = dd->adarrayghostedin[i]; 118 iarray_start = (char*)dd->adstartghostedin[i]; 119 itdof = dd->ghostedtdof; 120 dd->adarrayghostedin[i] = PETSC_NULL; 121 dd->adstartghostedin[i] = PETSC_NULL; 122 123 goto done; 124 } 125 } 126 xs = dd->Xs; 127 ys = dd->Ys; 128 zs = dd->Zs; 129 xm = dd->Xe-dd->Xs; 130 ym = dd->Ye-dd->Ys; 131 zm = dd->Ze-dd->Zs; 132 } else { 133 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 134 if (dd->adarrayin[i]) { 135 *iptr = dd->adarrayin[i]; 136 iarray_start = (char*)dd->adstartin[i]; 137 itdof = dd->tdof; 138 dd->adarrayin[i] = PETSC_NULL; 139 dd->adstartin[i] = PETSC_NULL; 140 141 goto done; 142 } 143 } 144 xs = dd->xs; 145 ys = dd->ys; 146 zs = dd->zs; 147 xm = dd->xe-dd->xs; 148 ym = dd->ye-dd->ys; 149 zm = dd->ze-dd->zs; 150 } 151 deriv_type_size = PetscADGetDerivTypeSize(); 152 153 switch (dd->dim) { 154 case 1: { 155 void *ptr; 156 itdof = xm; 157 158 ierr = PetscMalloc(xm*deriv_type_size,&iarray_start);CHKERRQ(ierr); 159 160 ptr = (void*)(iarray_start - xs*deriv_type_size); 161 *iptr = (void*)ptr; 162 break;} 163 case 2: { 164 void **ptr; 165 itdof = xm*ym; 166 167 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*deriv_type_size,&iarray_start);CHKERRQ(ierr); 168 169 ptr = (void**)(iarray_start + xm*ym*deriv_type_size - ys*sizeof(void*)); 170 for(j=ys;j<ys+ym;j++) { 171 ptr[j] = iarray_start + deriv_type_size*(xm*(j-ys) - xs); 172 } 173 *iptr = (void*)ptr; 174 break;} 175 case 3: { 176 void ***ptr,**bptr; 177 itdof = xm*ym*zm; 178 179 ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*deriv_type_size,&iarray_start);CHKERRQ(ierr); 180 181 ptr = (void***)(iarray_start + xm*ym*zm*deriv_type_size - zs*sizeof(void*)); 182 bptr = (void**)(iarray_start + xm*ym*zm*deriv_type_size + zm*sizeof(void**)); 183 184 for(i=zs;i<zs+zm;i++) { 185 ptr[i] = bptr + ((i-zs)*ym - ys); 186 } 187 for (i=zs; i<zs+zm; i++) { 188 for (j=ys; j<ys+ym; j++) { 189 ptr[i][j] = iarray_start + deriv_type_size*(xm*ym*(i-zs) + xm*(j-ys) - xs); 190 } 191 } 192 193 *iptr = (void*)ptr; 194 break;} 195 default: 196 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 197 } 198 199 done: 200 /* add arrays to the checked out list */ 201 if (ghosted) { 202 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 203 if (!dd->adarrayghostedout[i]) { 204 dd->adarrayghostedout[i] = *iptr ; 205 dd->adstartghostedout[i] = iarray_start; 206 dd->ghostedtdof = itdof; 207 break; 208 } 209 } 210 } else { 211 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 212 if (!dd->adarrayout[i]) { 213 dd->adarrayout[i] = *iptr ; 214 dd->adstartout[i] = iarray_start; 215 dd->tdof = itdof; 216 break; 217 } 218 } 219 } 220 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Too many DMDA ADIC arrays obtained"); 221 if (tdof) *tdof = itdof; 222 if (array_start) *(void**)array_start = iarray_start; 223 PetscFunctionReturn(0); 224 } 225 226 #undef __FUNCT__ 227 #define __FUNCT__ "DMDARestoreAdicArray" 228 /*@C 229 DMDARestoreAdicArray - Restores an array of derivative types for a DMDA 230 231 Input Parameter: 232 + da - information about my local patch 233 - ghosted - do you want arrays for the ghosted or nonghosted patch 234 235 Output Parameters: 236 + ptr - array data structured to be passed to ad_FormFunctionLocal() 237 . array_start - actual start of 1d array of all values that adiC can access directly 238 - tdof - total number of degrees of freedom represented in array_start 239 240 Level: advanced 241 242 .seealso: DMDAGetAdicArray() 243 244 @*/ 245 PetscErrorCode DMDARestoreAdicArray(DM da,PetscBool ghosted,void *ptr,void *array_start,PetscInt *tdof) 246 { 247 PetscInt i; 248 void **iptr = (void**)ptr,iarray_start = 0; 249 250 PetscFunctionBegin; 251 PetscValidHeaderSpecific(da,DM_CLASSID,1); 252 if (ghosted) { 253 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 254 if (dd->adarrayghostedout[i] == *iptr) { 255 iarray_start = dd->adstartghostedout[i]; 256 dd->adarrayghostedout[i] = PETSC_NULL; 257 dd->adstartghostedout[i] = PETSC_NULL; 258 break; 259 } 260 } 261 if (!iarray_start) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not find array in checkout list"); 262 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 263 if (!dd->adarrayghostedin[i]){ 264 dd->adarrayghostedin[i] = *iptr; 265 dd->adstartghostedin[i] = iarray_start; 266 break; 267 } 268 } 269 } else { 270 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 271 if (dd->adarrayout[i] == *iptr) { 272 iarray_start = dd->adstartout[i]; 273 dd->adarrayout[i] = PETSC_NULL; 274 dd->adstartout[i] = PETSC_NULL; 275 break; 276 } 277 } 278 if (!iarray_start) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not find array in checkout list"); 279 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 280 if (!dd->adarrayin[i]){ 281 dd->adarrayin[i] = *iptr; 282 dd->adstartin[i] = iarray_start; 283 break; 284 } 285 } 286 } 287 PetscFunctionReturn(0); 288 } 289 290 #undef __FUNCT__ 291 #define __FUNCT__ "ad_DAGetArray" 292 PetscErrorCode ad_DAGetArray(DM da,PetscBool ghosted,void *iptr) 293 { 294 PetscErrorCode ierr; 295 PetscFunctionBegin; 296 ierr = DMDAGetAdicArray(da,ghosted,iptr,0,0);CHKERRQ(ierr); 297 PetscFunctionReturn(0); 298 } 299 300 #undef __FUNCT__ 301 #define __FUNCT__ "ad_DARestoreArray" 302 PetscErrorCode ad_DARestoreArray(DM da,PetscBool ghosted,void *iptr) 303 { 304 PetscErrorCode ierr; 305 PetscFunctionBegin; 306 ierr = DMDARestoreAdicArray(da,ghosted,iptr,0,0);CHKERRQ(ierr); 307 PetscFunctionReturn(0); 308 } 309 310 #endif 311 312 #undef __FUNCT__ 313 #define __FUNCT__ "DMDAGetArray" 314 /*@C 315 DMDAGetArray - Gets a work array for a DMDA 316 317 Input Parameter: 318 + da - information about my local patch 319 - ghosted - do you want arrays for the ghosted or nonghosted patch 320 321 Output Parameters: 322 . vptr - array data structured 323 324 Note: The vector values are NOT initialized and may have garbage in them, so you may need 325 to zero them. 326 327 Level: advanced 328 329 .seealso: DMDARestoreArray(), DMDAGetAdicArray() 330 331 @*/ 332 PetscErrorCode DMDAGetArray(DM da,PetscBool ghosted,void *vptr) 333 { 334 PetscErrorCode ierr; 335 PetscInt j,i,xs,ys,xm,ym,zs,zm; 336 char *iarray_start; 337 void **iptr = (void**)vptr; 338 DM_DA *dd = (DM_DA*)da->data; 339 340 PetscFunctionBegin; 341 PetscValidHeaderSpecific(da,DM_CLASSID,1); 342 if (ghosted) { 343 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 344 if (dd->arrayghostedin[i]) { 345 *iptr = dd->arrayghostedin[i]; 346 iarray_start = (char*)dd->startghostedin[i]; 347 dd->arrayghostedin[i] = PETSC_NULL; 348 dd->startghostedin[i] = PETSC_NULL; 349 350 goto done; 351 } 352 } 353 xs = dd->Xs; 354 ys = dd->Ys; 355 zs = dd->Zs; 356 xm = dd->Xe-dd->Xs; 357 ym = dd->Ye-dd->Ys; 358 zm = dd->Ze-dd->Zs; 359 } else { 360 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 361 if (dd->arrayin[i]) { 362 *iptr = dd->arrayin[i]; 363 iarray_start = (char*)dd->startin[i]; 364 dd->arrayin[i] = PETSC_NULL; 365 dd->startin[i] = PETSC_NULL; 366 367 goto done; 368 } 369 } 370 xs = dd->xs; 371 ys = dd->ys; 372 zs = dd->zs; 373 xm = dd->xe-dd->xs; 374 ym = dd->ye-dd->ys; 375 zm = dd->ze-dd->zs; 376 } 377 378 switch (dd->dim) { 379 case 1: { 380 void *ptr; 381 382 ierr = PetscMalloc(xm*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 383 384 ptr = (void*)(iarray_start - xs*sizeof(PetscScalar)); 385 *iptr = (void*)ptr; 386 break;} 387 case 2: { 388 void **ptr; 389 390 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 391 392 ptr = (void**)(iarray_start + xm*ym*sizeof(PetscScalar) - ys*sizeof(void*)); 393 for(j=ys;j<ys+ym;j++) { 394 ptr[j] = iarray_start + sizeof(PetscScalar)*(xm*(j-ys) - xs); 395 } 396 *iptr = (void*)ptr; 397 break;} 398 case 3: { 399 void ***ptr,**bptr; 400 401 ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 402 403 ptr = (void***)(iarray_start + xm*ym*zm*sizeof(PetscScalar) - zs*sizeof(void*)); 404 bptr = (void**)(iarray_start + xm*ym*zm*sizeof(PetscScalar) + zm*sizeof(void**)); 405 for(i=zs;i<zs+zm;i++) { 406 ptr[i] = bptr + ((i-zs)*ym - ys); 407 } 408 for (i=zs; i<zs+zm; i++) { 409 for (j=ys; j<ys+ym; j++) { 410 ptr[i][j] = iarray_start + sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs); 411 } 412 } 413 414 *iptr = (void*)ptr; 415 break;} 416 default: 417 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 418 } 419 420 done: 421 /* add arrays to the checked out list */ 422 if (ghosted) { 423 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 424 if (!dd->arrayghostedout[i]) { 425 dd->arrayghostedout[i] = *iptr ; 426 dd->startghostedout[i] = iarray_start; 427 break; 428 } 429 } 430 } else { 431 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 432 if (!dd->arrayout[i]) { 433 dd->arrayout[i] = *iptr ; 434 dd->startout[i] = iarray_start; 435 break; 436 } 437 } 438 } 439 PetscFunctionReturn(0); 440 } 441 442 #undef __FUNCT__ 443 #define __FUNCT__ "DMDARestoreArray" 444 /*@C 445 DMDARestoreArray - Restores an array of derivative types for a DMDA 446 447 Input Parameter: 448 + da - information about my local patch 449 . ghosted - do you want arrays for the ghosted or nonghosted patch 450 - vptr - array data structured to be passed to ad_FormFunctionLocal() 451 452 Level: advanced 453 454 .seealso: DMDAGetArray(), DMDAGetAdicArray() 455 456 @*/ 457 PetscErrorCode DMDARestoreArray(DM da,PetscBool ghosted,void *vptr) 458 { 459 PetscInt i; 460 void **iptr = (void**)vptr,*iarray_start = 0; 461 DM_DA *dd = (DM_DA*)da->data; 462 463 PetscFunctionBegin; 464 PetscValidHeaderSpecific(da,DM_CLASSID,1); 465 if (ghosted) { 466 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 467 if (dd->arrayghostedout[i] == *iptr) { 468 iarray_start = dd->startghostedout[i]; 469 dd->arrayghostedout[i] = PETSC_NULL; 470 dd->startghostedout[i] = PETSC_NULL; 471 break; 472 } 473 } 474 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 475 if (!dd->arrayghostedin[i]){ 476 dd->arrayghostedin[i] = *iptr; 477 dd->startghostedin[i] = iarray_start; 478 break; 479 } 480 } 481 } else { 482 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 483 if (dd->arrayout[i] == *iptr) { 484 iarray_start = dd->startout[i]; 485 dd->arrayout[i] = PETSC_NULL; 486 dd->startout[i] = PETSC_NULL; 487 break; 488 } 489 } 490 for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) { 491 if (!dd->arrayin[i]){ 492 dd->arrayin[i] = *iptr; 493 dd->startin[i] = iarray_start; 494 break; 495 } 496 } 497 } 498 PetscFunctionReturn(0); 499 } 500 501 #undef __FUNCT__ 502 #define __FUNCT__ "DMDAGetAdicMFArray" 503 /*@C 504 DMDAGetAdicMFArray - Gets an array of derivative types for a DMDA for matrix-free ADIC. 505 506 Input Parameter: 507 + da - information about my local patch 508 - ghosted - do you want arrays for the ghosted or nonghosted patch? 509 510 Output Parameters: 511 + vptr - array data structured to be passed to ad_FormFunctionLocal() 512 . array_start - actual start of 1d array of all values that adiC can access directly (may be null) 513 - tdof - total number of degrees of freedom represented in array_start (may be null) 514 515 Notes: 516 The vector values are NOT initialized and may have garbage in them, so you may need 517 to zero them. 518 519 This routine returns the same type of object as the DMDAVecGetArray(), except its 520 elements are derivative types instead of PetscScalars. 521 522 Level: advanced 523 524 .seealso: DMDARestoreAdicMFArray(), DMDAGetArray(), DMDAGetAdicArray() 525 526 @*/ 527 PetscErrorCode DMDAGetAdicMFArray(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 528 { 529 PetscErrorCode ierr; 530 PetscInt j,i,xs,ys,xm,ym,zs,zm,itdof = 0; 531 char *iarray_start; 532 void **iptr = (void**)vptr; 533 DM_DA *dd = (DM_DA*)da->data; 534 535 PetscFunctionBegin; 536 PetscValidHeaderSpecific(da,DM_CLASSID,1); 537 if (ghosted) { 538 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 539 if (dd->admfarrayghostedin[i]) { 540 *iptr = dd->admfarrayghostedin[i]; 541 iarray_start = (char*)dd->admfstartghostedin[i]; 542 itdof = dd->ghostedtdof; 543 dd->admfarrayghostedin[i] = PETSC_NULL; 544 dd->admfstartghostedin[i] = PETSC_NULL; 545 546 goto done; 547 } 548 } 549 xs = dd->Xs; 550 ys = dd->Ys; 551 zs = dd->Zs; 552 xm = dd->Xe-dd->Xs; 553 ym = dd->Ye-dd->Ys; 554 zm = dd->Ze-dd->Zs; 555 } else { 556 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 557 if (dd->admfarrayin[i]) { 558 *iptr = dd->admfarrayin[i]; 559 iarray_start = (char*)dd->admfstartin[i]; 560 itdof = dd->tdof; 561 dd->admfarrayin[i] = PETSC_NULL; 562 dd->admfstartin[i] = PETSC_NULL; 563 564 goto done; 565 } 566 } 567 xs = dd->xs; 568 ys = dd->ys; 569 zs = dd->zs; 570 xm = dd->xe-dd->xs; 571 ym = dd->ye-dd->ys; 572 zm = dd->ze-dd->zs; 573 } 574 575 switch (dd->dim) { 576 case 1: { 577 void *ptr; 578 itdof = xm; 579 580 ierr = PetscMalloc(xm*2*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 581 582 ptr = (void*)(iarray_start - xs*2*sizeof(PetscScalar)); 583 *iptr = (void*)ptr; 584 break;} 585 case 2: { 586 void **ptr; 587 itdof = xm*ym; 588 589 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*2*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 590 591 ptr = (void**)(iarray_start + xm*ym*2*sizeof(PetscScalar) - ys*sizeof(void*)); 592 for(j=ys;j<ys+ym;j++) { 593 ptr[j] = iarray_start + 2*sizeof(PetscScalar)*(xm*(j-ys) - xs); 594 } 595 *iptr = (void*)ptr; 596 break;} 597 case 3: { 598 void ***ptr,**bptr; 599 itdof = xm*ym*zm; 600 601 ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*2*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 602 603 ptr = (void***)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) - zs*sizeof(void*)); 604 bptr = (void**)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) + zm*sizeof(void**)); 605 for(i=zs;i<zs+zm;i++) { 606 ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*); 607 } 608 for (i=zs; i<zs+zm; i++) { 609 for (j=ys; j<ys+ym; j++) { 610 ptr[i][j] = iarray_start + 2*sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs); 611 } 612 } 613 614 *iptr = (void*)ptr; 615 break;} 616 default: 617 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 618 } 619 620 done: 621 /* add arrays to the checked out list */ 622 if (ghosted) { 623 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 624 if (!dd->admfarrayghostedout[i]) { 625 dd->admfarrayghostedout[i] = *iptr ; 626 dd->admfstartghostedout[i] = iarray_start; 627 dd->ghostedtdof = itdof; 628 break; 629 } 630 } 631 } else { 632 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 633 if (!dd->admfarrayout[i]) { 634 dd->admfarrayout[i] = *iptr ; 635 dd->admfstartout[i] = iarray_start; 636 dd->tdof = itdof; 637 break; 638 } 639 } 640 } 641 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(((PetscObject)da)->comm,PETSC_ERR_ARG_WRONG,"Too many DMDA ADIC arrays obtained"); 642 if (tdof) *tdof = itdof; 643 if (array_start) *(void**)array_start = iarray_start; 644 PetscFunctionReturn(0); 645 } 646 647 #undef __FUNCT__ 648 #define __FUNCT__ "DMDAGetAdicMFArray4" 649 PetscErrorCode DMDAGetAdicMFArray4(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 650 { 651 PetscErrorCode ierr; 652 PetscInt j,i,xs,ys,xm,ym,itdof = 0; 653 char *iarray_start; 654 void **iptr = (void**)vptr; 655 DM_DA *dd = (DM_DA*)da->data; 656 657 PetscFunctionBegin; 658 PetscValidHeaderSpecific(da,DM_CLASSID,1); 659 if (ghosted) { 660 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 661 if (dd->admfarrayghostedin[i]) { 662 *iptr = dd->admfarrayghostedin[i]; 663 iarray_start = (char*)dd->admfstartghostedin[i]; 664 itdof = dd->ghostedtdof; 665 dd->admfarrayghostedin[i] = PETSC_NULL; 666 dd->admfstartghostedin[i] = PETSC_NULL; 667 668 goto done; 669 } 670 } 671 xs = dd->Xs; 672 ys = dd->Ys; 673 xm = dd->Xe-dd->Xs; 674 ym = dd->Ye-dd->Ys; 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 xm = dd->xe-dd->xs; 690 ym = dd->ye-dd->ys; 691 } 692 693 switch (dd->dim) { 694 case 2: { 695 void **ptr; 696 itdof = xm*ym; 697 698 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*5*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 699 700 ptr = (void**)(iarray_start + xm*ym*5*sizeof(PetscScalar) - ys*sizeof(void*)); 701 for(j=ys;j<ys+ym;j++) { 702 ptr[j] = iarray_start + 5*sizeof(PetscScalar)*(xm*(j-ys) - xs); 703 } 704 *iptr = (void*)ptr; 705 break;} 706 default: 707 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 708 } 709 710 done: 711 /* add arrays to the checked out list */ 712 if (ghosted) { 713 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 714 if (!dd->admfarrayghostedout[i]) { 715 dd->admfarrayghostedout[i] = *iptr ; 716 dd->admfstartghostedout[i] = iarray_start; 717 dd->ghostedtdof = itdof; 718 break; 719 } 720 } 721 } else { 722 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 723 if (!dd->admfarrayout[i]) { 724 dd->admfarrayout[i] = *iptr ; 725 dd->admfstartout[i] = iarray_start; 726 dd->tdof = itdof; 727 break; 728 } 729 } 730 } 731 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Too many DMDA ADIC arrays obtained"); 732 if (tdof) *tdof = itdof; 733 if (array_start) *(void**)array_start = iarray_start; 734 PetscFunctionReturn(0); 735 } 736 737 #undef __FUNCT__ 738 #define __FUNCT__ "DMDAGetAdicMFArray9" 739 PetscErrorCode DMDAGetAdicMFArray9(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 740 { 741 PetscErrorCode ierr; 742 PetscInt j,i,xs,ys,xm,ym,itdof = 0; 743 char *iarray_start; 744 void **iptr = (void**)vptr; 745 DM_DA *dd = (DM_DA*)da->data; 746 747 PetscFunctionBegin; 748 PetscValidHeaderSpecific(da,DM_CLASSID,1); 749 if (ghosted) { 750 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 751 if (dd->admfarrayghostedin[i]) { 752 *iptr = dd->admfarrayghostedin[i]; 753 iarray_start = (char*)dd->admfstartghostedin[i]; 754 itdof = dd->ghostedtdof; 755 dd->admfarrayghostedin[i] = PETSC_NULL; 756 dd->admfstartghostedin[i] = PETSC_NULL; 757 758 goto done; 759 } 760 } 761 xs = dd->Xs; 762 ys = dd->Ys; 763 xm = dd->Xe-dd->Xs; 764 ym = dd->Ye-dd->Ys; 765 } else { 766 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 767 if (dd->admfarrayin[i]) { 768 *iptr = dd->admfarrayin[i]; 769 iarray_start = (char*)dd->admfstartin[i]; 770 itdof = dd->tdof; 771 dd->admfarrayin[i] = PETSC_NULL; 772 dd->admfstartin[i] = PETSC_NULL; 773 774 goto done; 775 } 776 } 777 xs = dd->xs; 778 ys = dd->ys; 779 xm = dd->xe-dd->xs; 780 ym = dd->ye-dd->ys; 781 } 782 783 switch (dd->dim) { 784 case 2: { 785 void **ptr; 786 itdof = xm*ym; 787 788 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*10*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 789 790 ptr = (void**)(iarray_start + xm*ym*10*sizeof(PetscScalar) - ys*sizeof(void*)); 791 for(j=ys;j<ys+ym;j++) { 792 ptr[j] = iarray_start + 10*sizeof(PetscScalar)*(xm*(j-ys) - xs); 793 } 794 *iptr = (void*)ptr; 795 break;} 796 default: 797 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 798 } 799 800 done: 801 /* add arrays to the checked out list */ 802 if (ghosted) { 803 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 804 if (!dd->admfarrayghostedout[i]) { 805 dd->admfarrayghostedout[i] = *iptr ; 806 dd->admfstartghostedout[i] = iarray_start; 807 dd->ghostedtdof = itdof; 808 break; 809 } 810 } 811 } else { 812 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 813 if (!dd->admfarrayout[i]) { 814 dd->admfarrayout[i] = *iptr ; 815 dd->admfstartout[i] = iarray_start; 816 dd->tdof = itdof; 817 break; 818 } 819 } 820 } 821 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Too many DMDA ADIC arrays obtained"); 822 if (tdof) *tdof = itdof; 823 if (array_start) *(void**)array_start = iarray_start; 824 PetscFunctionReturn(0); 825 } 826 827 #undef __FUNCT__ 828 #define __FUNCT__ "DMDAGetAdicMFArrayb" 829 /*@C 830 DMDAGetAdicMFArrayb - Gets an array of derivative types for a DMDA for matrix-free ADIC. 831 832 Input Parameter: 833 + da - information about my local patch 834 - ghosted - do you want arrays for the ghosted or nonghosted patch? 835 836 Output Parameters: 837 + vptr - array data structured to be passed to ad_FormFunctionLocal() 838 . array_start - actual start of 1d array of all values that adiC can access directly (may be null) 839 - tdof - total number of degrees of freedom represented in array_start (may be null) 840 841 Notes: 842 The vector values are NOT initialized and may have garbage in them, so you may need 843 to zero them. 844 845 This routine returns the same type of object as the DMDAVecGetArray(), except its 846 elements are derivative types instead of PetscScalars. 847 848 Level: advanced 849 850 .seealso: DMDARestoreAdicMFArray(), DMDAGetArray(), DMDAGetAdicArray() 851 852 @*/ 853 PetscErrorCode DMDAGetAdicMFArrayb(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 854 { 855 PetscErrorCode ierr; 856 PetscInt j,i,xs,ys,xm,ym,zs,zm,itdof = 0; 857 char *iarray_start; 858 void **iptr = (void**)vptr; 859 DM_DA *dd = (DM_DA*)da->data; 860 PetscInt bs = dd->w,bs1 = bs+1; 861 862 PetscFunctionBegin; 863 PetscValidHeaderSpecific(da,DM_CLASSID,1); 864 if (ghosted) { 865 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 866 if (dd->admfarrayghostedin[i]) { 867 *iptr = dd->admfarrayghostedin[i]; 868 iarray_start = (char*)dd->admfstartghostedin[i]; 869 itdof = dd->ghostedtdof; 870 dd->admfarrayghostedin[i] = PETSC_NULL; 871 dd->admfstartghostedin[i] = PETSC_NULL; 872 873 goto done; 874 } 875 } 876 xs = dd->Xs; 877 ys = dd->Ys; 878 zs = dd->Zs; 879 xm = dd->Xe-dd->Xs; 880 ym = dd->Ye-dd->Ys; 881 zm = dd->Ze-dd->Zs; 882 } else { 883 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 884 if (dd->admfarrayin[i]) { 885 *iptr = dd->admfarrayin[i]; 886 iarray_start = (char*)dd->admfstartin[i]; 887 itdof = dd->tdof; 888 dd->admfarrayin[i] = PETSC_NULL; 889 dd->admfstartin[i] = PETSC_NULL; 890 891 goto done; 892 } 893 } 894 xs = dd->xs; 895 ys = dd->ys; 896 zs = dd->zs; 897 xm = dd->xe-dd->xs; 898 ym = dd->ye-dd->ys; 899 zm = dd->ze-dd->zs; 900 } 901 902 switch (dd->dim) { 903 case 1: { 904 void *ptr; 905 itdof = xm; 906 907 ierr = PetscMalloc(xm*bs1*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 908 909 ptr = (void*)(iarray_start - xs*bs1*sizeof(PetscScalar)); 910 *iptr = (void*)ptr; 911 break;} 912 case 2: { 913 void **ptr; 914 itdof = xm*ym; 915 916 ierr = PetscMalloc((ym+1)*sizeof(void*)+xm*ym*bs1*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 917 918 ptr = (void**)(iarray_start + xm*ym*bs1*sizeof(PetscScalar) - ys*sizeof(void*)); 919 for(j=ys;j<ys+ym;j++) { 920 ptr[j] = iarray_start + bs1*sizeof(PetscScalar)*(xm*(j-ys) - xs); 921 } 922 *iptr = (void*)ptr; 923 break;} 924 case 3: { 925 void ***ptr,**bptr; 926 itdof = xm*ym*zm; 927 928 ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*bs1*sizeof(PetscScalar),&iarray_start);CHKERRQ(ierr); 929 930 ptr = (void***)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) - zs*sizeof(void*)); 931 bptr = (void**)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) + zm*sizeof(void**)); 932 for(i=zs;i<zs+zm;i++) { 933 ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*); 934 } 935 for (i=zs; i<zs+zm; i++) { 936 for (j=ys; j<ys+ym; j++) { 937 ptr[i][j] = iarray_start + bs1*sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs); 938 } 939 } 940 941 *iptr = (void*)ptr; 942 break;} 943 default: 944 SETERRQ1(((PetscObject)da)->comm,PETSC_ERR_SUP,"Dimension %D not supported",dd->dim); 945 } 946 947 done: 948 /* add arrays to the checked out list */ 949 if (ghosted) { 950 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 951 if (!dd->admfarrayghostedout[i]) { 952 dd->admfarrayghostedout[i] = *iptr ; 953 dd->admfstartghostedout[i] = iarray_start; 954 dd->ghostedtdof = itdof; 955 break; 956 } 957 } 958 } else { 959 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 960 if (!dd->admfarrayout[i]) { 961 dd->admfarrayout[i] = *iptr ; 962 dd->admfstartout[i] = iarray_start; 963 dd->tdof = itdof; 964 break; 965 } 966 } 967 } 968 if (i == DMDA_MAX_AD_ARRAYS+1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Too many DMDA ADIC arrays obtained"); 969 if (tdof) *tdof = itdof; 970 if (array_start) *(void**)array_start = iarray_start; 971 PetscFunctionReturn(0); 972 } 973 974 #undef __FUNCT__ 975 #define __FUNCT__ "DMDARestoreAdicMFArray" 976 /*@C 977 DMDARestoreAdicMFArray - Restores an array of derivative types for a DMDA. 978 979 Input Parameter: 980 + da - information about my local patch 981 - ghosted - do you want arrays for the ghosted or nonghosted patch? 982 983 Output Parameters: 984 + ptr - array data structure to be passed to ad_FormFunctionLocal() 985 . array_start - actual start of 1d array of all values that adiC can access directly 986 - tdof - total number of degrees of freedom represented in array_start 987 988 Level: advanced 989 990 .seealso: DMDAGetAdicArray() 991 992 @*/ 993 PetscErrorCode DMDARestoreAdicMFArray(DM da,PetscBool ghosted,void *vptr,void *array_start,PetscInt *tdof) 994 { 995 PetscInt i; 996 void **iptr = (void**)vptr,*iarray_start = 0; 997 DM_DA *dd = (DM_DA*)da->data; 998 999 PetscFunctionBegin; 1000 PetscValidHeaderSpecific(da,DM_CLASSID,1); 1001 if (ghosted) { 1002 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1003 if (dd->admfarrayghostedout[i] == *iptr) { 1004 iarray_start = dd->admfstartghostedout[i]; 1005 dd->admfarrayghostedout[i] = PETSC_NULL; 1006 dd->admfstartghostedout[i] = PETSC_NULL; 1007 break; 1008 } 1009 } 1010 if (!iarray_start) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not find array in checkout list"); 1011 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1012 if (!dd->admfarrayghostedin[i]){ 1013 dd->admfarrayghostedin[i] = *iptr; 1014 dd->admfstartghostedin[i] = iarray_start; 1015 break; 1016 } 1017 } 1018 } else { 1019 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1020 if (dd->admfarrayout[i] == *iptr) { 1021 iarray_start = dd->admfstartout[i]; 1022 dd->admfarrayout[i] = PETSC_NULL; 1023 dd->admfstartout[i] = PETSC_NULL; 1024 break; 1025 } 1026 } 1027 if (!iarray_start) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not find array in checkout list"); 1028 for (i=0; i<DMDA_MAX_AD_ARRAYS; i++) { 1029 if (!dd->admfarrayin[i]){ 1030 dd->admfarrayin[i] = *iptr; 1031 dd->admfstartin[i] = iarray_start; 1032 break; 1033 } 1034 } 1035 } 1036 PetscFunctionReturn(0); 1037 } 1038 1039 #undef __FUNCT__ 1040 #define __FUNCT__ "admf_DAGetArray" 1041 PetscErrorCode admf_DAGetArray(DM da,PetscBool ghosted,void *iptr) 1042 { 1043 PetscErrorCode ierr; 1044 PetscFunctionBegin; 1045 ierr = DMDAGetAdicMFArray(da,ghosted,iptr,0,0);CHKERRQ(ierr); 1046 PetscFunctionReturn(0); 1047 } 1048 1049 #undef __FUNCT__ 1050 #define __FUNCT__ "admf_DARestoreArray" 1051 PetscErrorCode admf_DARestoreArray(DM da,PetscBool ghosted,void *iptr) 1052 { 1053 PetscErrorCode ierr; 1054 PetscFunctionBegin; 1055 ierr = DMDARestoreAdicMFArray(da,ghosted,iptr,0,0);CHKERRQ(ierr); 1056 PetscFunctionReturn(0); 1057 } 1058 1059