Skip to content

Commit

Permalink
Dynamic allocation of logging toggles (#126)
Browse files Browse the repository at this point in the history 
fixes #123

Co-authored-by: Pavel N. Krivitsky <p.krivitsky@unsw.edu.au>
  • Loading branch information
@AdrienLeGuillou @krivit
AdrienLeGuillou and krivit committed Oct 8, 2024
1 parent f84873b commit c24e7ac
Show file tree
Hide file tree
Showing 4 changed files with 61 additions and 25 deletions.
55 changes: 34 additions & 21 deletions src/MCMCDyn.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -52,20 +52,23 @@ SEXP MCMCDyn_wrapper(SEXP stateR, // ergm_state
memset(REAL(sample), 0, (asInteger(nsteps) + 1)*m->n_stats*sizeof(double));
memcpy(REAL(sample), s->stats, m->n_stats*sizeof(double));

SEXP difftime = PROTECT(allocVector(INTSXP, asInteger(log_changes) ? asInteger(maxchanges) : 0));
SEXP difftail = PROTECT(allocVector(INTSXP, asInteger(log_changes) ? asInteger(maxchanges) : 0));
SEXP diffhead = PROTECT(allocVector(INTSXP, asInteger(log_changes) ? asInteger(maxchanges) : 0));
SEXP diffto = PROTECT(allocVector(INTSXP, asInteger(log_changes) ? asInteger(maxchanges) : 0));
memset(INTEGER(difftime), 0, asInteger(log_changes) ? asInteger(maxchanges)*sizeof(int) : 0);
memset(INTEGER(difftail), 0, asInteger(log_changes) ? asInteger(maxchanges)*sizeof(int) : 0);
memset(INTEGER(diffhead), 0, asInteger(log_changes) ? asInteger(maxchanges)*sizeof(int) : 0);
memset(INTEGER(diffto), 0, asInteger(log_changes) ? asInteger(maxchanges)*sizeof(int) : 0);
kvint difftime, difftail, diffhead, diffto;
kv_init(difftime);
kv_init(difftail);
kv_init(diffhead);
kv_init(diffto);

// pre-allocate the 0th element for size
kv_push(int, difftime, 0);
kv_push(int, difftail, 0);
kv_push(int, diffhead, 0);
kv_push(int, diffto, 0);

SEXP status;
if(MHp) status = PROTECT(ScalarInteger(MCMCSampleDyn(s,
dur_inf,
REAL(eta),
asInteger(collect)?REAL(sample):NULL, asInteger(maxedges), asInteger(maxchanges), asInteger(log_changes), (Vertex *)INTEGER(difftime), (Vertex *)INTEGER(difftail), (Vertex *)INTEGER(diffhead), INTEGER(diffto),
asInteger(collect)?REAL(sample):NULL, asInteger(maxedges), asInteger(maxchanges), asInteger(log_changes), &difftime, &difftail, &diffhead, &diffto,
asInteger(nsteps), asInteger(min_MH_interval), asInteger(max_MH_interval), asReal(MH_pval), asReal(MH_interval_add), asInteger(burnin), asInteger(interval),
asInteger(verbose))));
else status = PROTECT(ScalarInteger(MCMCDyn_MH_FAILED));
Expand All @@ -81,10 +84,15 @@ SEXP MCMCDyn_wrapper(SEXP stateR, // ergm_state
SET_VECTOR_ELT(outl, 2, ErgmStateRSave(s));
}

SET_VECTOR_ELT(outl, 3, difftime);
SET_VECTOR_ELT(outl, 4, difftail);
SET_VECTOR_ELT(outl, 5, diffhead);
SET_VECTOR_ELT(outl, 6, diffto);
SET_VECTOR_ELT(outl, 3, PROTECT(kvint_to_SEXP(difftime)));
SET_VECTOR_ELT(outl, 4, PROTECT(kvint_to_SEXP(difftail)));
SET_VECTOR_ELT(outl, 5, PROTECT(kvint_to_SEXP(diffhead)));
SET_VECTOR_ELT(outl, 6, PROTECT(kvint_to_SEXP(diffto)));

kv_destroy(difftime);
kv_destroy(difftail);
kv_destroy(diffhead);
kv_destroy(diffto);

ErgmStateDestroy(s);
PutRNGstate(); /* Disable RNG before returning */
Expand All @@ -110,7 +118,7 @@ MCMCDynStatus MCMCSampleDyn(ErgmState *s,
int maxedges,
int maxchanges,
int log_changes,
Vertex *difftime, Vertex *difftail, Vertex *diffhead, int *diffto,
kvint *difftime, kvint *difftail, kvint *diffhead, kvint *diffto,
// MCMC settings.
unsigned int nsteps, unsigned int min_MH_interval, unsigned int max_MH_interval, double MH_pval, double MH_interval_add,
unsigned int burnin, unsigned int interval,
Expand Down Expand Up @@ -190,7 +198,12 @@ MCMCDynStatus MCMCSampleDyn(ErgmState *s,
}
}

if(log_changes) difftime[0]=difftail[0]=diffhead[0]=diffto[0]=nextdiffedge-1;
if(log_changes) {
kv_A(*difftime, 0) = nextdiffedge - 1;
kv_A(*difftail, 0) = nextdiffedge - 1;
kv_A(*diffhead, 0) = nextdiffedge - 1;
kv_A(*diffto, 0) = nextdiffedge - 1;
}
return MCMCDyn_OK;
}

Expand Down Expand Up @@ -223,7 +236,7 @@ MCMCDynStatus MCMCDyn1Step(ErgmState *s,
// Space for output.
double *stats,
unsigned int maxchanges, Edge *nextdiffedge,
Vertex *difftime, Vertex *difftail, Vertex *diffhead, int *diffto,
kvint *difftime, kvint *difftail, kvint *diffhead, kvint *diffto,
// MCMC settings.
unsigned int min_MH_interval, unsigned int max_MH_interval, double MH_pval, double MH_interval_add,
// Verbosity.
Expand Down Expand Up @@ -343,7 +356,7 @@ MCMCDynStatus MCMCDyn1Step_advance(ErgmState *s,
// Space for output.
double *stats,
unsigned int maxchanges, Edge *nextdiffedge,
Vertex *difftime, Vertex *difftail, Vertex *diffhead, int *diffto,
kvint *difftime, kvint *difftail, kvint *diffhead, kvint *diffto,
// Verbosity.
int verbose){
StoreDyadMapInt *discord = dur_inf->discord;
Expand All @@ -358,10 +371,10 @@ MCMCDynStatus MCMCDyn1Step_advance(ErgmState *s,
kh_foreach_key(discord, dyad,{
if(*nextdiffedge<maxchanges){
// and record the toggle.
if(difftime) difftime[*nextdiffedge] = t;
if(difftail) difftail[*nextdiffedge] = dyad.tail;
if(diffhead) diffhead[*nextdiffedge] = dyad.head;
if(diffto) diffto[*nextdiffedge] = GetEdge(dyad.tail,dyad.head,nwp);
if(difftime) kv_push(int, *difftime, t);
if(difftail) kv_push(int, *difftail, dyad.tail);
if(diffhead) kv_push(int, *diffhead, dyad.head);
if(diffto) kv_push(int, *diffto, GetEdge(dyad.tail,dyad.head,nwp));
(*nextdiffedge)++;
}else{
return(MCMCDyn_TOO_MANY_CHANGES);
Expand Down
9 changes: 5 additions & 4 deletions src/MCMCDyn.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@
#include "tergm_model.h"
#include "ergm_state.h"
#include "changestats_lasttoggle.h"
#include "kvint.h"

// TODO: This might be worth moving into a common "constants.h".
typedef enum MCMCDynStatus_enum {
Expand All @@ -34,10 +35,10 @@ MCMCDynStatus MCMCSampleDyn(ErgmState *s,
int maxedges,
int maxchanges,
int log_changes,
Vertex *difftime, Vertex *difftail, Vertex *diffhead, int *diffto,
kvint *difftime, kvint *difftail, kvint *diffhead, kvint *diffto,
// MCMC settings.
unsigned int nsteps, unsigned int min_MH_interval, unsigned int max_MH_interval, double MH_pval, double MH_interval_add,
unsigned int burnin, unsigned int interval,
unsigned int burnin, unsigned int interval,
// Verbosity.
int verbose);

Expand All @@ -47,7 +48,7 @@ MCMCDynStatus MCMCDyn1Step(ErgmState *s,
// Space for output.
double *stats,
unsigned int maxchanges, Edge *nextdiffedge,
Vertex *difftime, Vertex *difftail, Vertex *diffhead, int *diffto,
kvint *difftime, kvint *difftail, kvint *diffhead, kvint *diffto,
// MCMC settings.
unsigned int min_MH_interval, unsigned int max_MH_interval, double MH_pval, double MH_interval_add,
// Verbosity.
Expand All @@ -58,7 +59,7 @@ MCMCDynStatus MCMCDyn1Step_advance(ErgmState *s,
// Space for output.
double *stats,
unsigned int maxchanges, Edge *nextdiffedge,
Vertex *difftime, Vertex *difftail, Vertex *diffhead, int *diffto,
kvint *difftime, kvint *difftail, kvint *diffhead, kvint *diffto,
// Verbosity.
int verbose);
#endif
12 changes: 12 additions & 0 deletions src/kvint.c
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,12 @@
#include "kvint.h"

SEXP kvint_to_SEXP(kvint v) {
int size = kv_size(v);
SEXP out = PROTECT(allocVector(INTSXP, size));
int *out_ptr = INTEGER(out);
for (int i = 0; i < size; i++) {
out_ptr[i] = kv_A(v, i);
}
UNPROTECT(1);
return out;
}
10 changes: 10 additions & 0 deletions src/kvint.h
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,10 @@
#ifndef INCLUDE_SRC_KVINT_H_
#define INCLUDE_SRC_KVINT_H_

#include <Rinternals.h>
#include "ergm_kvec.h"

typedef kvec_t(int) kvint;
SEXP kvint_to_SEXP(kvint v);

#endif // INCLUDE_SRC_KVINT_H_

0 comments on commit c24e7ac

Please sign in to comment.