diff options
Diffstat (limited to 'src/lossdistrib.c')
| -rw-r--r-- | src/lossdistrib.c | 41 |
1 files changed, 41 insertions, 0 deletions
diff --git a/src/lossdistrib.c b/src/lossdistrib.c index b1b78c4..d8257bd 100644 --- a/src/lossdistrib.c +++ b/src/lossdistrib.c @@ -260,6 +260,47 @@ void lossdistrib_joint(const double *p, const double* pp, const int *np, free(qtemp); } +void joint_default_averagerecov_distrib(const double *p,const int *np, + const double *S, const int *N, + double *q) { + /* recursive algorithm with first order correction + computes jointly the loss and recovery distribution + p vector of default probabilities + np length of p + S vector of severities (should be same length as p) + N number of ticks in the grid + q the joint probability distribution */ + + double lu = 1./(*N-1); + q[0] = 1; + + double* newrecov = malloc(*N * sizeof(double)); + unsigned int* index = malloc(*N * sizeof(unsigned int)); + double* weights = malloc(*N * sizeof(double)); + q[0] = 1; + double temp; + for(int i = 0; i < (*np); i++) { + for(int j = i + 1; j > 0; j--) { + for(int n = 0; n < *N; n++) { + newrecov[n] = ( (j-1) * n * lu + (1-S[i])) / (double)j; + weights[n] = modf(newrecov[n] * (*N-1), &temp); + index[n] = (unsigned int) temp; + } + if(j < i+1) { + cblas_dscal(*N, 1-p[i], q+(*N)*j, 1); + } + for(int k = 0; k < *N; k++) { + q[j*(*N)+index[k]+1] += weights[k] * p[i] * q[(j-1)*(*N)+k]; + q[j*(*N)+index[k]] += (1-weights[k]) * p[i] * q[(j-1)*(*N)+k]; + } + } + cblas_dscal(*N, 1-p[i], q, 1); + } + free(newrecov); + free(index); +} + + void recovdist(const double *dp, const double *pp, const int *n, const double *w, const double *S, const int *N, double *q) { /* recursive algorithm with first order correction for computing |