version 1.3, 2018/03/30 04:43:16 |
version 1.16, 2019/10/21 05:37:20 |
|
|
/* $OpenXM: OpenXM/src/ox_gsl/ox_gsl.c,v 1.2 2018/03/29 11:52:18 takayama Exp $ |
/* $OpenXM: OpenXM/src/ox_gsl/ox_gsl.c,v 1.15 2018/06/08 00:03:43 takayama Exp $ |
*/ |
*/ |
|
|
#include <stdio.h> |
#include <stdio.h> |
|
|
#include <setjmp.h> |
#include <setjmp.h> |
#include <string.h> |
#include <string.h> |
#include <unistd.h> |
#include <unistd.h> |
|
#include <signal.h> |
#include <math.h> |
#include <math.h> |
#include "ox_gsl.h" |
#include "ox_gsl.h" |
#include "call_gsl.h" // need only when you bind call_gsl functions. |
#include "call_gsl.h" // need only when you bind call_gsl functions. |
|
#include "call_gsl_sf.h" |
|
|
OXFILE *fd_rw; |
OXFILE *fd_rw; |
|
|
Line 40 void gc_free(void *p,size_t size) |
|
Line 42 void gc_free(void *p,size_t size) |
|
void init_gc() |
void init_gc() |
{ GC_INIT(); |
{ GC_INIT(); |
mp_set_memory_functions(GC_malloc,gc_realloc,gc_free); |
mp_set_memory_functions(GC_malloc,gc_realloc,gc_free); |
|
init_dic(); // initialize ox_eval.c |
} |
} |
|
|
void initialize_stack() |
void initialize_stack() |
Line 110 int sm_mathcap() |
|
Line 113 int sm_mathcap() |
|
// CMO_DISTRIBUTED_POLYNOMIAL, |
// CMO_DISTRIBUTED_POLYNOMIAL, |
// CMO_RECURSIVE_POLYNOMIAL, |
// CMO_RECURSIVE_POLYNOMIAL, |
// CMO_POLYNOMIAL_IN_ONE_VARIABLE, |
// CMO_POLYNOMIAL_IN_ONE_VARIABLE, |
|
CMO_TREE, |
CMO_ERROR2, |
CMO_ERROR2, |
0}; |
0}; |
int available_sm_command[]={ |
int available_sm_command[]={ |
Line 154 cmo *make_error2(const char *reason,const char *fname, |
|
Line 158 cmo *make_error2(const char *reason,const char *fname, |
|
n = 5; |
n = 5; |
argv = (cmo **) GC_malloc(sizeof(cmo *)*n); |
argv = (cmo **) GC_malloc(sizeof(cmo *)*n); |
ms = (cmo *)new_cmo_string("Error"); argv[0] = ms; |
ms = (cmo *)new_cmo_string("Error"); argv[0] = ms; |
if (reason != NULL) s = (char *)GC_malloc(strlen(reason)+1); |
if (reason != NULL) {s = (char *)GC_malloc(strlen(reason)+1); strcpy(s,reason); |
else strcpy(s,""); |
}else strcpy(s,""); |
ms = (cmo *) new_cmo_string(s); argv[1] = ms; |
ms = (cmo *) new_cmo_string(s); argv[1] = ms; |
if (reason != NULL) s = (char *)GC_malloc(strlen(fname)+1); |
if (fname != NULL) {s = (char *)GC_malloc(strlen(fname)+1); strcpy(s,fname); |
else strcpy(s,""); |
}else strcpy(s,""); |
ms = (cmo *) new_cmo_string(s); argv[2] = ms; |
ms = (cmo *) new_cmo_string(s); argv[2] = ms; |
err = (cmo *)new_cmo_int32(line); argv[3] = err; |
err = (cmo *)new_cmo_int32(line); argv[3] = err; |
err = (cmo *)new_cmo_int32(code); argv[4] = err; |
err = (cmo *)new_cmo_int32(code); argv[4] = err; |
|
|
}else if (c->tag == CMO_ZERO) { |
}else if (c->tag == CMO_ZERO) { |
return(0); |
return(0); |
} |
} |
myhandler("get_i, not an integer",NULL,0,-1); |
myhandler("get_i: not an integer",NULL,0,-1); |
return 0; |
return 0; |
} |
} |
|
|
Line 208 double get_double() |
|
Line 212 double get_double() |
|
}else if (c->tag == CMO_ZZ) { |
}else if (c->tag == CMO_ZZ) { |
if ((mpz_cmp_si(mympz(c),(long int) 0x7fffffff)>0) || |
if ((mpz_cmp_si(mympz(c),(long int) 0x7fffffff)>0) || |
(mpz_cmp_si(mympz(c),(long int) -0x7fffffff)<0)) { |
(mpz_cmp_si(mympz(c),(long int) -0x7fffffff)<0)) { |
myhandler("get_double, out of int32",NULL,0,-1); |
myhandler("get_double: out of int32",NULL,0,-1); |
return(NAN); |
return(NAN); |
} |
} |
return( (double) mpz_get_si(((cmo_zz *)c)->mpz)); |
return( (double) mpz_get_si(((cmo_zz *)c)->mpz)); |
Line 217 double get_double() |
|
Line 221 double get_double() |
|
}else if (c->tag == CMO_ZERO) { |
}else if (c->tag == CMO_ZERO) { |
return(0); |
return(0); |
} |
} |
myhandler("get_double, not a double",NULL,0,-1); |
myhandler("get_double: not a double",NULL,0,-1); |
return(NAN); |
return(NAN); |
} |
} |
|
/* get_double() will be obsolted and will be replaced by cmo2double(c) */ |
|
double cmo2double(cmo *c) |
|
{ |
|
#define mympz(c) (((cmo_zz *)c)->mpz) |
|
if (c == NULL) c = pop(); |
|
if (c->tag == CMO_INT32) { |
|
return( (double) (((cmo_int32 *)c)->i) ); |
|
}else if (c->tag == CMO_IEEE_DOUBLE_FLOAT) { |
|
return (((cmo_double *)c)->d); // see ox_toolkit.h |
|
}else if (c->tag == CMO_ZZ) { |
|
if ((mpz_cmp_si(mympz(c),(long int) 0x7fffffff)>0) || |
|
(mpz_cmp_si(mympz(c),(long int) -0x7fffffff)<0)) { |
|
myhandler("get_double: out of int32",NULL,0,-1); |
|
return(NAN); |
|
} |
|
return( (double) mpz_get_si(((cmo_zz *)c)->mpz)); |
|
}else if (c->tag == CMO_NULL) { |
|
return(0); |
|
}else if (c->tag == CMO_ZERO) { |
|
return(0); |
|
} |
|
myhandler("cmo2double: not a double",NULL,0,-1); |
|
return(NAN); |
|
} |
|
|
void my_add_double() { |
void my_add_double() { |
double x,y; |
double x,y; |
Line 267 double *get_double_list(int *length) { |
|
Line 295 double *get_double_list(int *length) { |
|
} |
} |
return(d); |
return(d); |
} |
} |
|
/* get_double_list will be obsolted and will be replaced by cmo2double_list() */ |
|
double *cmo2double_list(int *length,cmo *c) { |
|
cmo *entry; |
|
cell *cellp; |
|
double *d; |
|
int n,i; |
|
if (c == NULL) c = pop(); |
|
if (c->tag != CMO_LIST) { |
|
// make_error2("get_double_list",NULL,0,-1); |
|
*length=-1; return(0); |
|
} |
|
n = *length = list_length((cmo_list *)c); |
|
d = (double *) GC_malloc(sizeof(double)*(*length+1)); |
|
cellp = list_first((cmo_list *)c); |
|
entry = cellp->cmo; |
|
for (i=0; i<n; i++) { |
|
if (Debug) { |
|
printf("entry[%d]=",i); print_cmo(entry); printf("\n"); |
|
} |
|
if (entry->tag == CMO_INT32) { |
|
d[i]=( (double) (((cmo_int32 *)entry)->i) ); |
|
}else if (entry->tag == CMO_IEEE_DOUBLE_FLOAT) { |
|
d[i]=((cmo_double *)entry)->d; |
|
}else if (entry->tag == CMO_ZZ) { |
|
d[i]=( (double) mpz_get_si(((cmo_zz *)entry)->mpz)); |
|
}else if (entry->tag == CMO_NULL) { |
|
d[i]= 0; |
|
}else { |
|
fprintf(stderr,"entries of the list should be int32 or zz or double\n"); |
|
*length = -1; |
|
myhandler("get_double_list",NULL,0,-1); |
|
return(NULL); |
|
} |
|
cellp = list_next(cellp); |
|
entry = cellp->cmo; |
|
} |
|
return(d); |
|
} |
void show_double_list() { |
void show_double_list() { |
int n; |
int n; |
double *d; |
double *d; |
Line 291 char *get_string() { |
|
Line 357 char *get_string() { |
|
return(NULL); |
return(NULL); |
} |
} |
|
|
|
void test_ox_eval() { |
|
cmo *c; |
|
double d=0; |
|
pop(); |
|
c=pop(); |
|
if (Debug) { |
|
ox_printf("cmo *c="); print_cmo(c); ox_printf("\n"); |
|
} |
|
init_dic(); |
|
register_entry("x",1.25); |
|
if (eval_cmo(c,&d) == 0) myhandler("eval_cmo failed",NULL,0,-1); |
|
push((cmo *)new_cmo_double(d)); |
|
} |
|
|
int sm_executeFunction() |
int sm_executeFunction() |
{ |
{ |
cmo_string *func = (cmo_string *)pop(); |
cmo_string *func = (cmo_string *)pop(); |
Line 298 int sm_executeFunction() |
|
Line 378 int sm_executeFunction() |
|
push(make_error2("sm_executeFunction, not CMO_STRING",NULL,0,-1)); |
push(make_error2("sm_executeFunction, not CMO_STRING",NULL,0,-1)); |
return -1; |
return -1; |
} |
} |
|
init_dic(); |
// Test functions |
// Test functions |
if (strcmp(func->s, "add_int32") == 0) { |
if (strcmp(func->s, "add_int32") == 0) { |
my_add_int32(); |
my_add_int32(); |
Line 307 int sm_executeFunction() |
|
Line 388 int sm_executeFunction() |
|
show_double_list(); |
show_double_list(); |
}else if (strcmp(func->s,"restart")==0) { |
}else if (strcmp(func->s,"restart")==0) { |
pop(); restart(); |
pop(); restart(); |
|
}else if (strcmp(func->s,"test_ox_eval")==0) { |
|
test_ox_eval(); |
// The following functions are defined in call_gsl.c |
// The following functions are defined in call_gsl.c |
}else if (strcmp(func->s,"gsl_sf_lngamma_complex_e")==0) { |
}else if (strcmp(func->s,"gsl_sf_lngamma_complex_e")==0) { |
call_gsl_sf_lngamma_complex_e(); |
call_gsl_sf_lngamma_complex_e(); |
|
}else if (strcmp(func->s,"gsl_integration_qags")==0) { |
|
call_gsl_integration_qags(); |
|
}else if (strcmp(func->s,"gsl_monte_plain_integrate")==0) { |
|
call_gsl_monte_plain_miser_vegas_integrate(0); |
|
}else if (strcmp(func->s,"gsl_monte_miser_integrate")==0) { |
|
call_gsl_monte_plain_miser_vegas_integrate(1); |
|
}else if (strcmp(func->s,"gsl_monte_vegas_integrate")==0) { |
|
call_gsl_monte_plain_miser_vegas_integrate(2); |
|
}else if (strcmp(func->s,"gsl_odeiv_step_rk4")==0) { |
|
call_gsl_odeiv_step("rk4"); |
|
}else if (strcmp(func->s,"gsl_sf_gamma_inc")==0) { |
|
call_gsl_sf_gamma_inc(); |
}else { |
}else { |
push(make_error2("sm_executeFunction, unknown function",NULL,0,-1)); |
push(make_error2("sm_executeFunction, unknown function",NULL,0,-1)); |
return -1; |
return -1; |
Line 382 void myhandler(const char *reason,const char *file,int |
|
Line 477 void myhandler(const char *reason,const char *file,int |
|
fprintf(fp,"%d\n",line); |
fprintf(fp,"%d\n",line); |
if (file != NULL) fprintf(fp,"%s\n",file); else fprintf(fp,"file?\n"); |
if (file != NULL) fprintf(fp,"%s\n",file); else fprintf(fp,"file?\n"); |
if (reason != NULL) fprintf(fp,"%s\n",reason); else fprintf(fp,"reason?\n"); |
if (reason != NULL) fprintf(fp,"%s\n",reason); else fprintf(fp,"reason?\n"); |
fflush(NULL); fclose(fp); // BUG. the contents are deleted after it is closed. |
fflush(NULL); fclose(fp); |
// m = make_error2(reason,file,line,gsl_errno); |
// m = make_error2(reason,file,line,gsl_errno); |
// send_ox_cmo(fd_rw, m); ox_flush(fd_rw); |
// send_ox_cmo(fd_rw, m); ox_flush(fd_rw); |
// send error packet even it is not asked. Todo, OK? --> no |
// send error packet even it is not asked. Todo, OK? --> no |
restart(); |
restart(); |
} |
} |
|
|
void push_error_from_file() { |
void push_error_from_file() { |
FILE *fp; |
FILE *fp; |
#define BUF_SIZE 1024 |
#define BUF_SIZE 1024 |
Line 400 void push_error_from_file() { |
|
Line 494 void push_error_from_file() { |
|
cmo *m; |
cmo *m; |
fprintf(stderr,"push_error_from_file()\n"); |
fprintf(stderr,"push_error_from_file()\n"); |
sprintf(logname,"/tmp/ox_gsl-%d.txt",(int) getpid()); |
sprintf(logname,"/tmp/ox_gsl-%d.txt",(int) getpid()); |
fp = fopen(logname,"w"); |
fp = fopen(logname,"r"); |
if (fp == NULL) return; |
if (fp == NULL) { |
|
fprintf(stderr,"open %s is failed\n",logname); return; |
|
} |
fgets(cmd,BUF_SIZE-2,fp); sscanf(cmd,"%d",&gsl_errno); |
fgets(cmd,BUF_SIZE-2,fp); sscanf(cmd,"%d",&gsl_errno); |
fgets(cmd,BUF_SIZE-2,fp); sscanf(cmd,"%d",&line); |
fgets(cmd,BUF_SIZE-2,fp); sscanf(cmd,"%d",&line); |
fgets(file,BUF_SIZE-2,fp); |
#define remove_newline(s) {char *tmp_pos; if ((tmp_pos=strchr(s,'\n')) != NULL) *tmp_pos = '\0';} |
fgets(reason,BUF_SIZE-2,fp); |
fgets(file,BUF_SIZE-2,fp); remove_newline(file); |
|
fgets(reason,BUF_SIZE-2,fp); remove_newline(reason); |
fclose(fp); |
fclose(fp); |
m = make_error2(reason,file,line,gsl_errno); |
m = make_error2(reason,file,line,gsl_errno); |
push(m); |
push(m); |
|
|
fd_rw = oxf_open(3); |
fd_rw = oxf_open(3); |
oxf_determine_byteorder_server(fd_rw); |
oxf_determine_byteorder_server(fd_rw); |
} |
} |
|
#if defined(__CYGWIN__) |
|
void *mysignal(int sig,void (*handler)(int m)); |
|
mysignal(SIGUSR1,usr1_handler); |
|
#else |
signal(SIGUSR1,usr1_handler); |
signal(SIGUSR1,usr1_handler); |
|
#endif |
|
|
while(1) { |
while(1) { |
receive(); |
receive(); |
} |
} |
return(0); |
return(0); |
|
} |
|
|
|
cmo *element_of_at(cmo *list,int k) { |
|
int length; |
|
static cmo * saved_list = NULL; |
|
static cmo **dic; |
|
int i; |
|
cell *cellp; |
|
if (list == NULL) { |
|
ox_printf("element_of_at: list is NULL.\n"); |
|
return( (cmo *)NULL); |
|
} |
|
if (list->tag != CMO_LIST) { |
|
ox_printf("element_of_at: list is not list.\n"); |
|
return((cmo *)NULL); |
|
} |
|
length = list_length((cmo_list *)list); |
|
if ((k < 0) || (k >= length)) { |
|
ox_printf("element_of_at: out of bound length=%d, k=%d.\n",length,k); |
|
return((cmo *)NULL); |
|
} |
|
if (list == saved_list) return(dic[k]); |
|
saved_list = list; |
|
dic = (cmo **)GC_malloc(sizeof(cmo *)*(length+1)); |
|
if (dic == NULL) return((cmo *)NULL); // no more memory. |
|
cellp = list_first((cmo_list *)list); |
|
for (i=0; i<length; i++) { |
|
dic[i] = cellp->cmo; |
|
cellp = list_next(cellp); |
|
} |
|
return(dic[k]); |
|
} |
|
|
|
int get_length(cmo *c) { |
|
if (c->tag != CMO_LIST) { |
|
return(-1); |
|
} |
|
return(list_length((cmo_list *)c)); |
} |
} |