[BACK]Return to imach.c CVS log [TXT] [DIR] Up to [Local Repository] / imach / src

Diff for /imach/src/imach.c between versions 1.361 and 1.362

version 1.361, 2024/05/12 20:29:32 version 1.362, 2024/06/28 08:00:31
Line 1 Line 1
 /* $Id$  /* $Id$
   $State$    $State$
   $Log$    $Log$
     Revision 1.362  2024/06/28 08:00:31  brouard
     Summary: 0.99s6
   
     * imach.c (Module): s6 errors with age*age (harmless).
   
   Revision 1.361  2024/05/12 20:29:32  brouard    Revision 1.361  2024/05/12 20:29:32  brouard
   Summary: Version 0.99s5    Summary: Version 0.99s5
   
Line 1531  char *endptr; Line 1536  char *endptr;
 long lval;  long lval;
 double dval;  double dval;
   
   /* This for praxis gegen */
     /* int prin=1; */
     double h0=0.25;
     double macheps;
     double ffmin;
   
 #define NR_END 1  #define NR_END 1
 #define FREE_ARG char*  #define FREE_ARG char*
 #define FTOL 1.0e-10  #define FTOL 1.0e-10
Line 3231  L1:  /* L1 or try loop */ Line 3242  L1:  /* L1 or try loop */
       if (k > 0) *d2 = 0;        if (k > 0) *d2 = 0;
    }     }
 #ifdef DEBUGPRAX  #ifdef DEBUGPRAX
    printf(" bebe end of min x1=%14.8e fx=%14.8e d2=%14.8e\n",*x1, fx, *d2);     printf(" bebe end of min x1 might be very wrong x1=%14.8e fx=%14.8e d2=%14.8e\n",*x1, fx, *d2);
 #endif  #endif
    if (*d2 <= small_windows) *d2 = small_windows;     if (*d2 <= small_windows) *d2 = small_windows;
    *x1 = x2; fx = fm;     *x1 = x2; fx = fm;
Line 3707  mloop: Line 3718  mloop:
    printf("praxis4 macheps=%14g h=%14g step=%14g small=%14g t=%14g\n",macheps,h, h0,small_windows, t);      printf("praxis4 macheps=%14g h=%14g step=%14g small=%14g t=%14g\n",macheps,h, h0,small_windows, t); 
 #endif  #endif
    /* min(0, 2, &d[0], &s, fx, 0); /\* mac heps not global *\/ */     /* min(0, 2, &d[0], &s, fx, 0); /\* mac heps not global *\/ */
    minny(1, 2, &d[1], &s, fx, 0); /* mac heps not global */     minny(1, 2, &d[1], &s, fx, 0); /* mac heps not global it seems that fx doesn't correspond to f(s=*x1) */
 #ifdef DEBUGPRAX  #ifdef DEBUGPRAX
    printf("praxis5 macheps=%14g h=%14g looks at sign of s=%14g fx=%14g\n",macheps,h, s,fx);      printf("praxis5 macheps=%14g h=%14g looks at sign of s=%14g fx=%14g\n",macheps,h, s,fx); 
 #endif  #endif
Line 4218  void powell(double p[], double **xi, int Line 4229  void powell(double p[], double **xi, int
         printf("  + age*age  ");          printf("  + age*age  ");
         fprintf(ficlog,"  + age*age  ");          fprintf(ficlog,"  + age*age  ");
     }      }
     for(j=1;j <=ncovmodel-2;j++){      for(j=1;j <=ncovmodel-2-nagesqr;j++){
       if(Typevar[j]==0) {        if(Typevar[j]==0) {
         printf("  +      V%d  ",Tvar[j]);          printf("  +      V%d  ",Tvar[j]);
         fprintf(ficlog,"  +      V%d  ",Tvar[j]);          fprintf(ficlog,"  +      V%d  ",Tvar[j]);
Line 6593  void mlikeli(FILE *ficres,double p[], in Line 6604  void mlikeli(FILE *ficres,double p[], in
 #else  /* FLATSUP */  #else  /* FLATSUP */
 /*  powell(p,xi,npar,ftol,&iter,&fret,func);*/  /*  powell(p,xi,npar,ftol,&iter,&fret,func);*/
 /*   praxis ( t0, h0, n, prin, x, beale_f ); */  /*   praxis ( t0, h0, n, prin, x, beale_f ); */
   int prin=1;    /* int prin=1; */
   double h0=0.25;    /* double h0=0.25; */
   double macheps;    /* double macheps; */
   double fmin;    /* double fmin; */
   macheps=pow(16.0,-13.0);    macheps=pow(16.0,-13.0);
 /* #include "praxis.h" */  /* #include "praxis.h" */
   /* Be careful that praxis start at x[0] and powell start at p[1] */    /* Be careful that praxis start at x[0] and powell start at p[1] */
Line 6606  printf("Praxis Gegenfurtner \n"); Line 6617  printf("Praxis Gegenfurtner \n");
 fprintf(ficlog, "Praxis  Gegenfurtner\n");fflush(ficlog);  fprintf(ficlog, "Praxis  Gegenfurtner\n");fflush(ficlog);
 /* praxis ( ftol, h0, npar, prin, p1, func ); */  /* praxis ( ftol, h0, npar, prin, p1, func ); */
   /* fmin = praxis(1.e-5,macheps, h, n, prin, x, func); */    /* fmin = praxis(1.e-5,macheps, h, n, prin, x, func); */
   fmin = praxis(ftol,macheps, h0, npar, prin, p, func);    ffmin = praxis(ftol,macheps, h0, npar, prin, p, func);
 printf("End Praxis\n");  printf("End Praxis\n");
 #endif  /* FLATSUP */  #endif  /* FLATSUP */
   
Line 8602  void  concatwav(int wav[], int **dh, int Line 8613  void  concatwav(int wav[], int **dh, int
    double ***gradg, ***trgradg; /**< for var eij */     double ***gradg, ***trgradg; /**< for var eij */
    double **gradgp, **trgradgp; /**< for var p point j */     double **gradgp, **trgradgp; /**< for var p point j */
    double *gpp, *gmp; /**< for var p point j */     double *gpp, *gmp; /**< for var p point j */
    double **varppt; /**< for var e.. nlstate+1 to nlstate+ndeath */     double **varppt; /**< for var p.3 p.death nlstate+1 to nlstate+ndeath */
    double ***p3mat;     double ***p3mat;
    double age,agelim, hf;     double age,agelim, hf;
    /* double ***mobaverage; */     /* double ***mobaverage; */
Line 12283  double gompertz(double x[]) Line 12294  double gompertz(double x[])
          A=-x[1]/(x[2])*(exp(x[2]*(agecens[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp)));           A=-x[1]/(x[2])*(exp(x[2]*(agecens[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp)));
        } else if (cens[i] == 0){         } else if (cens[i] == 0){
         A=-x[1]/(x[2])*(exp(x[2]*(agedc[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp)))          A=-x[1]/(x[2])*(exp(x[2]*(agedc[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp)))
           +log(x[1]/YEARM) +x[2]*(agedc[i]-agegomp)+log(YEARM);            +log(fabs(x[1])/YEARM) +x[2]*(agedc[i]-agegomp)+log(YEARM);
           /* +log(x[1]/YEARM) +x[2]*(agedc[i]-agegomp)+log(YEARM); */  /* To be seen */
       } else        } else
          printf("Gompertz cens[%d] neither 1 nor 0\n",i);           printf("Gompertz cens[%d] neither 1 nor 0\n",i);
       /*if (wav[i] > 1 && agecens[i] > 15) {*/ /* ??? */        /*if (wav[i] > 1 && agecens[i] > 15) {*/ /* ??? */
Line 15843  Interval (in months) between two waves: Line 15855  Interval (in months) between two waves:
   flatdir=ivector(1,npar);     flatdir=ivector(1,npar); 
   for (j=1;j<=npar;j++) flatdir[j]=0;     for (j=1;j<=npar;j++) flatdir[j]=0; 
 #endif /*LINMINORIGINAL */  #endif /*LINMINORIGINAL */
      powell(p,ximort,NDIM,ftol,&iter,&fret,gompertz);      /* powell(p,ximort,NDIM,ftol,&iter,&fret,gompertz); */
 #endif      /* double h0=0.25; */
     macheps=pow(16.0,-13.0);
     printf("Praxis Gegenfurtner mle=%d\n",mle);
     fprintf(ficlog, "Praxis  Gegenfurtner mle=%d\n", mle);fflush(ficlog);
      /* ffmin = praxis(ftol,macheps, h0, npar, prin, p, gompertz); */
     /* For the Gompertz we use only two parameters */
     int _npar=2;
      ffmin = praxis(ftol,macheps, h0, _npar, 4, p, gompertz);
     printf("End Praxis\n");
     fclose(ficrespow);      fclose(ficrespow);
 #ifdef LINMINORIGINAL  #ifdef LINMINORIGINAL
 #else  #else
Line 15984  Please run with mle=-1 to get a correct Line 16004  Please run with mle=-1 to get a correct
       fprintf(ficlog,"  + age*age  ");        fprintf(ficlog,"  + age*age  ");
       fprintf(fichtm, "<th>+ age*age</th>");        fprintf(fichtm, "<th>+ age*age</th>");
     }      }
     for(j=1;j <=ncovmodel-2;j++){      for(j=1;j <=ncovmodel-2-nagesqr;j++){
       if(Typevar[j]==0) {        if(Typevar[j]==0) {
         printf("  +      V%d  ",Tvar[j]);          printf("  +      V%d  ",Tvar[j]);
         fprintf(ficres,"  +      V%d  ",Tvar[j]);          fprintf(ficres,"  +      V%d  ",Tvar[j]);
Line 16055  Please run with mle=-1 to get a correct Line 16075  Please run with mle=-1 to get a correct
         fprintf(ficlog,"  + age*age  ");          fprintf(ficlog,"  + age*age  ");
         fprintf(fichtm, "<th>+ age*age</th>");          fprintf(fichtm, "<th>+ age*age</th>");
       }        }
       for(j=1;j <=ncovmodel-2;j++){        for(j=1;j <=ncovmodel-2-nagesqr;j++){
         if(Typevar[j]==0) {          if(Typevar[j]==0) {
           printf("  +      V%d  ",Tvar[j]);            printf("  +      V%d  ",Tvar[j]);
           fprintf(fichtm, "<th>+ V%d</th>",Tvar[j]);            fprintf(fichtm, "<th>+ V%d</th>",Tvar[j]);
Removed from v.1.361  
changed lines
  Added in v.1.362

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>