1
$\begingroup$

The documentation of PETSc shows an example of usage:

  PetscRandomCreate(PETSC_COMM_SELF,&r);
  PetscRandomSetType(r,PETSCRAND48);
  PetscRandomGetValue(r,&value1);
  PetscRandomGetValueReal(r,&value2);
  PetscRandomDestroy(&r);

In the example random 1, we can find also the use of PetscRandomSeed.

However, it seems not to be enough for attaching a PetscRandom variable to an object.

The following code works well for any number of procs:

#undef __FUNCT__
#define __FUNCT__ "TSInjectCreate"
PetscErrorCode  TSInjectCreate(MPI_Comm comm,TSInject *ininject){
  PetscErrorCode ierr;
  size_t         t;
  PetscInt       i;
  PetscReal       r;
  TSInject       inject;

  PetscFunctionBegin;
  PetscValidPointer(ininject,1);
  *ininject = NULL;
  ierr                       = PetscClassIdRegister("TSInject",&TSINJECT_CLASSID);CHKERRQ(ierr);
  ierr = PetscHeaderCreate(inject,TSINJECT_CLASSID,"TSInject","Injection of SDC","TS",comm,TSInjectDestroy, TSInjectView);CHKERRQ(ierr);
  ierr = PetscRandomCreate( comm,&inject->ran);CHKERRQ(ierr);
  ierr = PetscRandomSetType(inject->ran, PETSCRAND48);CHKERRQ(ierr);
  ierr = PetscRandomSeed(inject->ran);CHKERRQ(ierr);
  ierr = PetscRandomGetValue(inject->ran, &r);CHKERRQ(ierr);
  *ininject = inject;

  PetscFunctionReturn(0);
}

But as soon as I use PetscRandomGetValue(inject->ran, &r); in another function, it works only when mpi is launched with one proc.

Why?

EDIT:

The error comes exactly (and weirdly) from the comparison:

PetscErrorCode  TSInjection(TSInject inject, Vec X)
{
  PetscReal       r,err;
  Vec             Y;
  PetscErrorCode ierr;

  PetscFunctionBegin;
  ierr = PetscRandomGetValue(inject->ran, &r);CHKERRQ(ierr);
ierr = VecView(X, 0);CHKERRQ(ierr);    // works
  if(r  > inject->proba){ PetscFunctionReturn(0);}
        ierr = VecView(X, 0);CHKERRQ(ierr);      // stucks the execution
 ...          
  PetscFunctionReturn(0);
}

If I replace r by 0.1 for example, everything works fine. What is all the more surprising is that displaying r doesn't stuck the execution.

Improve this question
$\endgroup$
3
  • $\begingroup$ Is there a reason your inject object should survive passing out of scope? $\endgroup$
    – origimbo
    Commented Mar 22, 2016 at 23:17
  • $\begingroup$ Yes the reason is that it receives parameters and these parameters will choose which function should be called. In this initialization function, non important lines were deleted. $\endgroup$
    – guhur
    Commented Mar 23, 2016 at 13:43
  • $\begingroup$ Moreover, the problem occurs on several machines. Linux/mac, same issue. $\endgroup$
    – guhur
    Commented Mar 23, 2016 at 15:32

1 Answer 1

Reset to default
1
$\begingroup$

To solve this issue, a work-around is too employ a DIY random function. In my context, I have given up with PetscRandom to write:

PetscErrorCode TSInjectRandInt(TSInject inject, PetscInt *r)  // RAND_MAX assumed to be 32767
{
    PetscFunctionBegin;
    inject->next = inject->next * 1103515245 + 12345;
    *r = PetscAbs((PetscInt)(inject->next/65536) % 32768);
    PetscFunctionReturn(0);
}


PetscErrorCode TSInjectRandReal(TSInject inject, PetscReal *r)  {
  PetscInt i;

  PetscFunctionBegin;
  TSInjectRandInt(inject, &i);
  *r = ((PetscReal)i)/32767.0;
  PetscFunctionReturn(0);
}

PetscErrorCode TSInjectRandSeed(TSInject inject, PetscInt seed)
{
    PetscFunctionBegin;
    inject->next = seed;
    PetscFunctionReturn(0);
}

However, one should notice that RAND_MAX is fixed and that the generator is clearly very easy. It has the only advantage of working.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.