← Back to branch summary

~ubuntu-branches/ubuntu/vivid/lme4/vivid

« back to all changes in this revision

Viewing changes to src/predModule.cpp

  • Committer: Package Import Robot
  • Author(s): Dirk Eddelbuettel
  • Date: 2014-08-28 14:04:58 UTC
  • mfrom: (1.2.25)
  • Revision ID: package-import@ubuntu.com-20140828140458-88n08mg8ec84w4w8
Tags: 1.1-7-1
* New upstream release

* Upload was delayed for seven weeks it took for the new Build-Depends
  r-cran-nlopr to get into unstable.  NEW queue is usually much quicker,

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/predModule.cpp
57
57
        setTheta(d_theta);          // starting values into Lambda
58
58
        d_L.cholmod().final_ll = 1; // force an LL' decomposition
59
59
        updateLamtUt();
60
 
        d_L.analyzePattern(d_LamtUt); // perform symbolic analysis
 
60
        d_L.analyzePattern(d_LamtUt * d_LamtUt.transpose()); // perform symbolic analysis
61
61
        if (d_L.info() != Eigen::Success)
62
62
            throw runtime_error("CholeskyDecomposition.analyzePattern failed");
63
63
    }
249
249
    
250
250
    // using a point so as to detect NULL
251
251
    void merPredD::updateDecomp(const MatrixXd* xPenalty) {  // update L, RZX and RX
 
252
        int debug=0;
 
253
        
 
254
        if (debug) Rcpp::Rcout << "start updateDecomp" << std::endl;
252
255
        updateL();
 
256
        if (debug) {
 
257
            Rcpp::Rcout << "updateDecomp 2: " << 
 
258
                d_RZX.cols() << " " << d_RZX.rows() << " " <<
 
259
                d_LamtUt.cols() << " " << d_LamtUt.rows() << " " <<
 
260
                d_V.cols() << " " << d_V.rows() << " " <<
 
261
                std::endl;
 
262
        }
 
263
        if (d_LamtUt.cols() != d_V.rows()) {
 
264
            ::Rf_warning("dimension mismatch in updateDecomp()");
 
265
            // Rcpp::Rcout << "WARNING: dimension mismatch in updateDecomp(): " <<
 
266
            // " LamtUt=" << d_LamtUt.rows() << "x" << d_LamtUt.cols() << 
 
267
            // "; V=" << d_V.rows() << "x" << d_V.cols() << " " <<
 
268
            // std::endl;
 
269
        }
253
270
        d_RZX         = d_LamtUt * d_V;
 
271
        if (debug) Rcpp::Rcout << "updateDecomp 3" << std::endl;
254
272
        if (d_p > 0) {
255
273
            d_L.solveInPlace(d_RZX, CHOLMOD_P);
256
274
            d_L.solveInPlace(d_RZX, CHOLMOD_L);
257
 
            
 
275
            if (debug) Rcpp::Rcout << "updateDecomp 4" << std::endl;
258
276
            MatrixXd      VtVdown(d_VtV);
259
277
            
260
278
            if (xPenalty == NULL)
262
280
            else {
263
281
                d_RX.compute(VtVdown.selfadjointView<Eigen::Upper>().rankUpdate(d_RZX.adjoint(), -1).rankUpdate(*xPenalty, 1));
264
282
            }
 
283
            if (debug) Rcpp::Rcout << "updateDecomp 5" << std::endl;
265
284
            if (d_RX.info() != Eigen::Success)
266
285
                ::Rf_error("Downdated VtV is not positive definite");
267
286
            d_ldRX2         = 2. * d_RX.matrixLLT().diagonal().array().abs().log().sum();
 
287
            if (debug) Rcpp::Rcout << "updateDecomp 6" << std::endl;
268
288
        }
269
289
    }
270
290