1d9262e54SJed Brown 2c6db04a5SJed Brown #include <petscsys.h> /*I "petscsys.h" I*/ 3d9262e54SJed Brown 48bf1f09cSShri Abhyankar #if defined(PETSC_USE_DEBUG) 5d9262e54SJed Brown 61f46d60fSShri Abhyankar #if defined(PETSC_HAVE_PTHREADCLASSES) 7997ce2baSShri Abhyankar #if defined(PETSC_PTHREAD_LOCAL) 8997ce2baSShri Abhyankar PETSC_PTHREAD_LOCAL PetscStack *petscstack = 0; 91f46d60fSShri Abhyankar #endif 10fe89fe5aSShri Abhyankar #else 117087cfbeSBarry Smith PetscStack *petscstack = 0; 12fe89fe5aSShri Abhyankar #endif 13d9262e54SJed Brown 1415681b3cSBarry Smith 1515681b3cSBarry Smith #if defined(PETSC_HAVE_AMS) 1615681b3cSBarry Smith #include <petscviewerams.h> 1715681b3cSBarry Smith 18*ec957eceSBarry Smith static AMS_Memory amsmemstack = NULL; 1915681b3cSBarry Smith 20d9262e54SJed Brown #undef __FUNCT__ 2115681b3cSBarry Smith #define __FUNCT__ "PetscStackAMSGrantAccess" 2215681b3cSBarry Smith /*@C 2315681b3cSBarry Smith PetscStackAMSGrantAccess - Grants access of the PETSc stack frames to the AMS publisher 2415681b3cSBarry Smith 2515681b3cSBarry Smith Collective on PETSC_COMM_WORLD? 2615681b3cSBarry Smith 2715681b3cSBarry Smith Level: developer 2815681b3cSBarry Smith 2915681b3cSBarry Smith Concepts: publishing object 3015681b3cSBarry Smith 3115681b3cSBarry Smith Developers Note: Cannot use PetscFunctionBegin/Return() or PetscStackCallAMS() since it may be used within those routines 3215681b3cSBarry Smith 3315681b3cSBarry Smith .seealso: PetscObjectSetName(), PetscObjectAMSViewOff(), PetscObjectAMSTakeAccess() 3415681b3cSBarry Smith 3515681b3cSBarry Smith @*/ 3615681b3cSBarry Smith void PetscStackAMSGrantAccess(void) 37d9262e54SJed Brown { 38*ec957eceSBarry Smith if (amsmemstack) { 39*ec957eceSBarry Smith AMS_Unlock_Memory(amsmemstack); 4015681b3cSBarry Smith } 41d9262e54SJed Brown } 42d9262e54SJed Brown 43d9262e54SJed Brown #undef __FUNCT__ 4415681b3cSBarry Smith #define __FUNCT__ "PetscStackAMSTakeAccess" 4515681b3cSBarry Smith /*@C 4615681b3cSBarry Smith PetscStackAMSTakeAccess - Takes access of the PETSc stack frames to the AMS publisher 4715681b3cSBarry Smith 4815681b3cSBarry Smith Collective on PETSC_COMM_WORLD? 4915681b3cSBarry Smith 5015681b3cSBarry Smith Level: developer 5115681b3cSBarry Smith 5215681b3cSBarry Smith Concepts: publishing object 5315681b3cSBarry Smith 5415681b3cSBarry Smith Developers Note: Cannot use PetscFunctionBegin/Return() or PetscStackCallAMS() since it may be used within those routines 5515681b3cSBarry Smith 5615681b3cSBarry Smith .seealso: PetscObjectSetName(), PetscObjectAMSViewOff(), PetscObjectAMSTakeAccess() 5715681b3cSBarry Smith 5815681b3cSBarry Smith @*/ 5915681b3cSBarry Smith void PetscStackAMSTakeAccess(void) 60d9262e54SJed Brown { 61*ec957eceSBarry Smith if (amsmemstack) { 62*ec957eceSBarry Smith AMS_Lock_Memory(amsmemstack); 6315681b3cSBarry Smith } 6415681b3cSBarry Smith } 6515681b3cSBarry Smith 6615681b3cSBarry Smith PetscErrorCode PetscStackViewAMS(void) 6715681b3cSBarry Smith { 6815681b3cSBarry Smith AMS_Memory mem; 6915681b3cSBarry Smith PetscStack* petscstackp; 7015681b3cSBarry Smith 7115681b3cSBarry Smith petscstackp = (PetscStack*)PetscThreadLocalGetValue(petscstack); 72*ec957eceSBarry Smith PetscStackCallAMS(AMS_Memory_Create,("Stack",&mem)); 73*ec957eceSBarry Smith PetscStackCallAMS(AMS_New_Field,(mem,"functions",petscstackp->function,10,AMS_READ,AMS_STRING)); 74*ec957eceSBarry Smith PetscStackCallAMS(AMS_New_Field,(mem,"current size",&petscstackp->currentsize,1,AMS_READ,AMS_INT)); 7515681b3cSBarry Smith amsmemstack = mem; 7615681b3cSBarry Smith return 0; 7715681b3cSBarry Smith } 7815681b3cSBarry Smith 7915681b3cSBarry Smith #undef __FUNCT__ 8015681b3cSBarry Smith #define __FUNCT__ "PetscStackAMSViewOff" 8115681b3cSBarry Smith PetscErrorCode PetscStackAMSViewOff(void) 8215681b3cSBarry Smith { 8315681b3cSBarry Smith PetscErrorCode ierr; 8415681b3cSBarry Smith 85d9262e54SJed Brown PetscFunctionBegin; 86*ec957eceSBarry Smith if (!amsmemstack) PetscFunctionReturn(0); 87*ec957eceSBarry Smith ierr = AMS_Memory_Destroy(amsmemstack);CHKERRQ(ierr); 88*ec957eceSBarry Smith amsmemstack = NULL; 89d9262e54SJed Brown PetscFunctionReturn(0); 90d9262e54SJed Brown } 91d9262e54SJed Brown 9215681b3cSBarry Smith #endif 9315681b3cSBarry Smith 947d5f7e0cSShri Abhyankar PetscErrorCode PetscStackCreate(void) 9574b43855SShri Abhyankar { 9674b43855SShri Abhyankar PetscStack *petscstack_in; 97dbf62e16SBarry Smith if (PetscStackActive()) return 0; 9874b43855SShri Abhyankar 9974b43855SShri Abhyankar petscstack_in = (PetscStack*)malloc(sizeof(PetscStack)); 10074b43855SShri Abhyankar petscstack_in->currentsize = 0; 101047240e1SBarry Smith PetscThreadLocalSetValue((PetscThreadKey*)&petscstack,petscstack_in); 10215681b3cSBarry Smith 10315681b3cSBarry Smith #if defined(PETSC_HAVE_AMS) 10415681b3cSBarry Smith { 10515681b3cSBarry Smith PetscBool flg = PETSC_FALSE; 10615681b3cSBarry Smith PetscOptionsHasName(NULL,"-stack_view",&flg); 10715681b3cSBarry Smith if (flg) PetscStackViewAMS(); 10815681b3cSBarry Smith } 10915681b3cSBarry Smith #endif 11074b43855SShri Abhyankar return 0; 11174b43855SShri Abhyankar } 11274b43855SShri Abhyankar 11315681b3cSBarry Smith 114d9262e54SJed Brown #undef __FUNCT__ 115d9262e54SJed Brown #define __FUNCT__ "PetscStackView" 116639ff905SBarry Smith PetscErrorCode PetscStackView(FILE *file) 117d9262e54SJed Brown { 118d9262e54SJed Brown int i; 1191f46d60fSShri Abhyankar PetscStack *petscstackp; 120d9262e54SJed Brown 1211f46d60fSShri Abhyankar petscstackp = (PetscStack*)PetscThreadLocalGetValue(petscstack); 122639ff905SBarry Smith if (!file) file = PETSC_STDOUT; 123d9262e54SJed Brown 124d9262e54SJed Brown if (file == PETSC_STDOUT) { 125d9262e54SJed Brown (*PetscErrorPrintf)("Note: The EXACT line numbers in the stack are not available,\n"); 126d9262e54SJed Brown (*PetscErrorPrintf)(" INSTEAD the line number of the start of the function\n"); 127d9262e54SJed Brown (*PetscErrorPrintf)(" is given.\n"); 128a297a907SKarl Rupp for (i=petscstackp->currentsize-1; i>=0; i--) (*PetscErrorPrintf)("[%d] %s line %d %s%s\n",PetscGlobalRank,petscstackp->function[i],petscstackp->line[i],petscstackp->directory[i],petscstackp->file[i]); 129d9262e54SJed Brown } else { 130d9262e54SJed Brown fprintf(file,"Note: The EXACT line numbers in the stack are not available,\n"); 131d9262e54SJed Brown fprintf(file," INSTEAD the line number of the start of the function\n"); 132d9262e54SJed Brown fprintf(file," is given.\n"); 133a297a907SKarl Rupp for (i=petscstackp->currentsize-1; i>=0; i--) fprintf(file,"[%d] %s line %d %s%s\n",PetscGlobalRank,petscstackp->function[i],petscstackp->line[i],petscstackp->directory[i],petscstackp->file[i]); 134d9262e54SJed Brown } 135d9262e54SJed Brown return 0; 136d9262e54SJed Brown } 137d9262e54SJed Brown 1387d5f7e0cSShri Abhyankar PetscErrorCode PetscStackDestroy(void) 13974b43855SShri Abhyankar { 140dbf62e16SBarry Smith if (PetscStackActive()) { 1411f46d60fSShri Abhyankar PetscStack *petscstack_in; 1421f46d60fSShri Abhyankar petscstack_in = (PetscStack*)PetscThreadLocalGetValue(petscstack); 14374b43855SShri Abhyankar free(petscstack_in); 144047240e1SBarry Smith PetscThreadLocalSetValue((PetscThreadKey*)&petscstack,(PetscStack*)0); 14561d886c9SShri Abhyankar PetscThreadLocalDestroy(petscstack); /* Deletes pthread_key if it was used */ 1467d5f7e0cSShri Abhyankar } 147d9262e54SJed Brown return 0; 148d9262e54SJed Brown } 149d9262e54SJed Brown 150d9262e54SJed Brown #undef __FUNCT__ 151d9262e54SJed Brown #define __FUNCT__ "PetscStackCopy" 152d9262e54SJed Brown /* PetscFunctionBegin; so that make rule checkbadPetscFunctionBegin works */ 1537087cfbeSBarry Smith PetscErrorCode PetscStackCopy(PetscStack *sint,PetscStack *sout) 154d9262e54SJed Brown { 155d9262e54SJed Brown int i; 156d9262e54SJed Brown 157a297a907SKarl Rupp if (!sint) sout->currentsize = 0; 158a297a907SKarl Rupp else { 159d9262e54SJed Brown for (i=0; i<sint->currentsize; i++) { 160d9262e54SJed Brown sout->function[i] = sint->function[i]; 161d9262e54SJed Brown sout->file[i] = sint->file[i]; 162d9262e54SJed Brown sout->directory[i] = sint->directory[i]; 163d9262e54SJed Brown sout->line[i] = sint->line[i]; 164a8d2bbe5SBarry Smith sout->petscroutine[i] = sint->petscroutine[i]; 165d9262e54SJed Brown } 166d9262e54SJed Brown sout->currentsize = sint->currentsize; 167d9262e54SJed Brown } 168d9262e54SJed Brown return 0; 169d9262e54SJed Brown } 170d9262e54SJed Brown 171d9262e54SJed Brown #undef __FUNCT__ 172d9262e54SJed Brown #define __FUNCT__ "PetscStackPrint" 173d9262e54SJed Brown /* PetscFunctionBegin; so that make rule checkbadPetscFunctionBegin works */ 1747087cfbeSBarry Smith PetscErrorCode PetscStackPrint(PetscStack *sint,FILE *fp) 175d9262e54SJed Brown { 176d9262e54SJed Brown int i; 177d9262e54SJed Brown 178d9262e54SJed Brown if (!sint) return(0); 179a297a907SKarl Rupp for (i=sint->currentsize-2; i>=0; i--) fprintf(fp," [%d] %s() line %d in %s%s\n",PetscGlobalRank,sint->function[i],sint->line[i],sint->directory[i],sint->file[i]); 180d9262e54SJed Brown return 0; 181d9262e54SJed Brown } 182d9262e54SJed Brown 183d9262e54SJed Brown #else 18415681b3cSBarry Smith 185d9262e54SJed Brown #undef __FUNCT__ 186d9262e54SJed Brown #define __FUNCT__ "PetscStackCreate" 1877087cfbeSBarry Smith PetscErrorCode PetscStackCreate(void) 188d9262e54SJed Brown { 189d9262e54SJed Brown PetscFunctionBegin; 190d9262e54SJed Brown PetscFunctionReturn(0); 191d9262e54SJed Brown } 192d9262e54SJed Brown #undef __FUNCT__ 193d9262e54SJed Brown #define __FUNCT__ "PetscStackView" 194639ff905SBarry Smith PetscErrorCode PetscStackView(FILE *file) 195d9262e54SJed Brown { 196d9262e54SJed Brown PetscFunctionBegin; 197d9262e54SJed Brown PetscFunctionReturn(0); 198d9262e54SJed Brown } 199d9262e54SJed Brown #undef __FUNCT__ 200d9262e54SJed Brown #define __FUNCT__ "PetscStackDestroy" 2017087cfbeSBarry Smith PetscErrorCode PetscStackDestroy(void) 202d9262e54SJed Brown { 203d9262e54SJed Brown PetscFunctionBegin; 204d9262e54SJed Brown PetscFunctionReturn(0); 205d9262e54SJed Brown } 206d9262e54SJed Brown 2076eb7abf7SJed Brown #if defined(PETSC_HAVE_AMS) /* AMS stack functions do nothing in optimized mode */ 2086eb7abf7SJed Brown void PetscStackAMSGrantAccess(void) {} 2096eb7abf7SJed Brown void PetscStackAMSTakeAccess(void) {} 2106eb7abf7SJed Brown 2116eb7abf7SJed Brown PetscErrorCode PetscStackViewAMS(void) 2126eb7abf7SJed Brown { 2136eb7abf7SJed Brown return 0; 2146eb7abf7SJed Brown } 2156eb7abf7SJed Brown 2166eb7abf7SJed Brown #undef __FUNCT__ 2176eb7abf7SJed Brown #define __FUNCT__ "PetscStackAMSViewOff" 2186eb7abf7SJed Brown PetscErrorCode PetscStackAMSViewOff(void) 2196eb7abf7SJed Brown { 2206eb7abf7SJed Brown PetscFunctionBegin; 2216eb7abf7SJed Brown PetscFunctionReturn(0); 2226eb7abf7SJed Brown } 2236eb7abf7SJed Brown #endif 2246eb7abf7SJed Brown 225d9262e54SJed Brown #endif 226d9262e54SJed Brown 227