Resources

elegant threading 2

0

Jason,

As requested, here is the updated, much messier version of my code.  Fixing the seeds did not keep it from crashing.  Something else is going on.  Ugh.

Thanks,

Chuck

Here are the call key variable defs and the thread statements:

thrd3_dif = thrd_dif*3;
Zfull1 = zeros(thrd3_dif[1],2); Zfull2 = zeros(thrd3_dif[2],2); Zfull3 = zeros(thrd3_dif[3],2);
Zbi1 = zeros(thrd3_dif[1],1); Zbi2 = zeros(thrd3_dif[2],1); Zbi3 = zeros(thrd3_dif[3],1);
Cap_xi1 = zeros(thrd3_dif[1],1); Cap_xi2 = zeros(thrd3_dif[2],1); Cap_xi3 = zeros(thrd3_dif[3],1);
bi1 = zeros(rc_dif[1],1); bi2 = zeros(rc_dif[2],1); bi3 = zeros(rc_dif[3],1);
OD1_drw = zeros(thrd_dif[1],1); OD2_drw = zeros(thrd_dif[2],1); OD3_drw = zeros(thrd_dif[3],1);
DAP1_drw = zeros(thrd_dif[1],1); DAP2_drw = zeros(thrd_dif[2],1); DAP3_drw = zeros(thrd_dif[3],1);

for g (1,11000,1);
//print "g = " g;

 @ draw updated values of the mean parameters @
 threadStat {cXiVX1,cXiVw1,Zfull1,Zbi1,Cap_xi1,bi1,OD1_drw,DAP1_drw,p0,d1,pr} = reg_parms(thrd_idx[1],(thrd_idx[2]-1),rc_dif[1],seed[1],p0, d0, d1, pr);
 threadStat {cXiVX2,cXiVw2,Zfull2,Zbi2,Cap_xi2,bi2,OD2_drw,DAP2_drw,p0,d1,pr} = reg_parms(thrd_idx[2],(thrd_idx[3]-1),rc_dif[2],seed[2],p0, d0, d1, pr);
 threadStat {cXiVX3,cXiVw3,Zfull3,Zbi3,Cap_xi3,bi3,OD3_drw,DAP3_drw,p0,d1,pr} = reg_parms(thrd_idx[3],(thrd_idx[4]-1),rc_dif[3],seed[3],p0, d0, d1, pr);
 threadJoin;

cXiVX = cXiVX1+cXiVX2+cXiVX3;
cXiVw = cXiVw1+cXiVw2+cXiVw3;
d1.Zfull = Zfull1|Zfull2|Zfull3;
d1.Zbi = Zbi1|Zbi2|Zbi3;
d1.Cap_xi = Cap_xi1|Cap_xi2|Cap_xi3;
p0.bi_sp1 = bi1|bi2|bi3;
pr.OD_drw = OD1_drw|OD2_drw|OD3_drw;
pr.DAP_drw = DAP1_drw|DAP2_drw|DAP3_drw;

 and her is the proc.  There are no global writes to highlight this time  :) :

proc(11) = reg_parms(_strt,_end,_rc_dif,_seed,struct parameter p0, struct data0 d0, struct data1 d1, struct prob pr);
 local _cXiVX,_cXiVw,imz,_imz,imz2,_imzj,imzj,i,j,ij,_ij,ci,ci3,i1,i2,D_i,cZiVZ,cZiVw,Zmat,ZiVZ,ZiVw,icZiVZ,cicZiVZ,bi_hat;
 local bnd,l_bnd,rnd,Xmat,XiVX,XiVw,U01,_rw,_rw3,_Zfull,_Zbi,_Cap_xi,_bi_sp1,_ci,_OD_drw,_DAP_drw,bi_std;
 
 _cXiVX = zeros(XODcol,XODcol);
 _cXiVw = zeros(XODcol,1);
 
 _rw = (_end - _strt + 1);
 _rw3 = 3*_rw;
 _Zfull = zeros(_rw3,2);
 _Zbi = zeros(_rw3,1);
 _Cap_xi = zeros(_rw3,1);
 _bi_sp1 = zeros(_rc_dif,2);
 _OD_drw = zeros(_rw,1);
 _DAP_drw = zeros(_rw,1);
 bnd = zeros(2,1);
 l_bnd = zeros(2,1);
 rnd = zeros(2,1);
 
 _imz=0; imz=3*_strt-3;
 _ci=0;
 imz2=2*_strt-2;
 for i (_strt,_end,1);
  imz=imz+3;
  _imz=_imz+3;
  
  @ start calculations for bi_sp1[z,.] @
  if cust_chg[i] == 1;
   ci = cust_id[i];
   _ci=_ci+1;
   ci3=3*ci;
   imz2=imz2+2;
   _imzj=_imz-3; imzj=imz-3;
   D_i = d1.Dfull[imz2-1:imz2,.];
   cZiVZ = zeros(2,2);
   cZiVw = zeros(2,1);
   for j (0,c_cnt[ci]-1,1);
     imzj=imzj+3; _imzj=_imzj+3;
    ij = i+j;
    _ij = ij-_strt+1;
    if d0.DAP_only[ij] == 1;
     {i1} = offer(pr.CumOP_O);
     _OD_drw[_ij] = pr.iOD_fee[i1];
     _DAP_drw[_ij] = dap_fees[ij];
     Zmat = ((pr.iOD_fee[i1]/d1.deposits[ij])~0)|(0~d1.dap_dep[ij])|(pr.iOD_fee[i1]/(d1.deposits[ij])~d1.dap_dep[ij]);
     _Zfull[_imzj-2:_imzj,.] = Zmat;
     
    elseif d0.OD_only[ij] == 1;
     {i1} = offer(pr.CumOP_D);
     _OD_drw[_ij] = od_fees[ij];
     _DAP_drw[_ij] = pr.iDAP_fee[i1];
     Zmat = (d1.od_dep[ij]~0)|(0~(pr.iDAP_fee[i1]/d1.deposits[ij]))|(d1.od_dep[ij]~(pr.iDAP_fee[i1]/d1.deposits[ij]));
     _Zfull[_imzj-2:_imzj,.] = Zmat;
     
    elseif d0.ODDAP[ij] == 1;
     _OD_drw[_ij] = od_fees[ij];
     _DAP_drw[_ij] = dap_fees[ij];
     Zmat = d1.Zfull[imzj-2:imzj,.];
     _Zfull[_imzj-2:_imzj,.] = Zmat;
    elseif d0.DAP_only[ij] == 0 and d0.OD_only[ij] == 0 and d0.ODDAP[ij] == 0;
     {i1} = offer(pr.CumOP_O);
     {i2} = offer(pr.CumOP_D);
     _OD_drw[_ij] = pr.iOD_fee[i1];
     _DAP_drw[_ij] = pr.iDAP_fee[i2];
     Zmat = ((pr.iOD_fee[i1]/d1.deposits[ij])~0)|(0~(pr.iDAP_fee[i2]/d1.deposits[ij]))|
      ((pr.iOD_fee[i1]/d1.deposits[ij])~(pr.iDAP_fee[i2]/d1.deposits[ij]));
     _Zfull[_imzj-2:_imzj,.] = Zmat;
     endif;
    ZiVZ = Zmat'p0.iV_sp1*Zmat;
    ZiVw = Zmat'p0.iV_sp1*(d0.w_sp1[imzj-2:imzj] - d1.Xfull[imzj-2:imzj,.]*p0.b_sp1 - p0.Cap_xi[ci3-2:ci3]);
    cZiVZ = cZiVZ + ZiVZ;
    cZiVw = cZiVw + ZiVw;
    if j == c_cnt[ci]-1;
     icZiVZ = invpd(cZiVZ + p0.iH_sp1);
     bi_hat = icZiVZ*(cZiVw + p0.iH_sp1*D_i*p0.thet_sp1);
     
     cicZiVZ = chol(icZiVZ)';
     {bi_std,_seed} = rndlcn(2,1,_seed);
     @ bi_sp1[ci,.] = (bi_hat + cicZiVZ*bi_std)';  without the GHK @
     
     @ GHK Simulator @
     {U01,_seed} = rndlcu(2,1,_seed);
     bnd[1] = -bi_hat[1]/cicZiVZ[1,1];
     l_bnd[1] = (-5 - bi_hat[1])/cicZiVZ[1,1];
     rnd[1] = cdfni(U01[1]*(cdfn(bnd[1]) - cdfn(l_bnd[1])) + cdfn(l_bnd[1]));
     bnd[2] = (-bi_hat[2] - cicZiVZ[2,1]*rnd[1])/cicZiVZ[2,2];
     l_bnd[2] = (-5 - bi_hat[2] - cicZiVZ[2,1]*rnd[1])/cicZiVZ[2,2];
     rnd[2] = cdfni(U01[2]*(cdfn(bnd[2]) - cdfn(l_bnd[2])) + cdfn(l_bnd[2]));
     _bi_sp1[_ci,.] = (bi_hat + cicZiVZ*rnd)';
   endif;
   endfor; @ j loop @
  endif;
  _Zbi[_imz-2:_imz] = _Zfull[_imz-2:_imz,.]*_bi_sp1[_ci,.]';
  _Cap_xi[_imz-2:_imz] = p0.Cap_xi[ci3-2:ci3];
 
  @ end calculations for bi_sp1[ci,.] @

  Xmat = d1.Xfull[imz-2:imz,.];
  XiVX = Xmat'p0.iV_sp1*Xmat;
  XiVw = Xmat'p0.iV_sp1*(d0.w_sp1[imz-2:imz] - _Zbi[_imz-2:_imz] - p0.Cap_xi[ci3-2:ci3]);
  _cXiVX = _cXiVX + XiVX;
  _cXiVw = _cXiVw + XiVw;
 endfor; @ i loop @
 retp(_cXiVX,_cXiVw,_Zfull,_Zbi,_Cap_xi,_bi_sp1,_OD_drw,_DAP_drw,p0,d1,pr);
endp;

 

asked September 11, 2013

1 Answer

0

Hello Chuck,

Each of the threadStat statements below is assigning to p0, d1 and pr. These are global assignments happening inside of the threads.

//draw updated values of the mean parameters
threadStat { cXiVX1,cXiVw1,Zfull1,Zbi1,Cap_xi1,bi1,
             OD1_drw,DAP1_drw,p0,d1,pr} = reg_parms(thrd_idx[1],(thrd_idx[2]-1),rc_dif[1],seed[1],p0, d0, d1, pr); 
threadStat {cXiVX2,cXiVw2,Zfull2,Zbi2,Cap_xi2,bi2,
            OD2_drw,DAP2_drw,p0,d1,pr} = reg_parms(thrd_idx[2],(thrd_idx[3]-1),rc_dif[2],seed[2],p0, d0, d1, pr); 
threadStat {cXiVX3,cXiVw3,Zfull3,Zbi3,Cap_xi3,bi3,
            OD3_drw,DAP3_drw,p0,d1,pr} = reg_parms(thrd_idx[3],(thrd_idx[4]-1),rc_dif[3],seed[3],p0, d0, d1, pr); 
threadJoin;