source: trunk/ADOL-C/examples/additional_examples/openmp_exam/liborpar.cpp @ 708

Last change on this file since 708 was 708, checked in by kulshres, 3 years ago

Merge branch 'master' of 'gitclone' into 'svn'

The following changes have been merged:

commit e2291bde44a282a133894b0db350aeb0b92a87db
Author: Mladen Banovic <mladenbanovic2705@…>
Date: Fri Jul 8 10:15:51 2016 +0200

Add methods getNumLiveVar and getNumDir in adtl.h, change counter type in FOR_I_EQ_0_LT_NUMDIR macro to size_t (instead of int). Update chunk size of BOOST pool in adouble_tl.cpp according to adouble::numDir.

commit 2ffb294465b973bfd4bf1f73d84478f8233c0d2f
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Thu Jun 23 12:32:14 2016 +0200

implement missing ref_eq_mult_p und ref_eq_min_p in ho_rev.c

somehow these were left out when parameters were being implemented.

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

commit 8cf0e5c1bd36f1dcf3be72cd67de631b2e1d0ee6
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Thu Jun 23 12:31:04 2016 +0200

make sure the result is the last locint written in trace for each operation

since we're trying to generate ascii traces in the future, we'll need this
convention that the last location is the result, and previous locations
are arguments. This has been the case for almost all operations anyway
except for a few new one's that I wrote without keeping this in mind.

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

commit 9ae0ff220f37463f2ed85cafc8a626c24e472f2f
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Tue Jun 21 14:16:27 2016 +0200

on some compilers newer boost interferes with AC_FUNC_MALLOC test

so do AC_FUNC_MALLOC and AC_FUNC_REALLOC as usual and check for boost
library later.

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

commit b746f620772cc8cce53e8f350adc6281279caf72
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Mon Jun 20 15:32:22 2016 +0200

make Klaus Röbenack's name UTF-8 instead of ISO-8859-1

These are the only places where we're not simple ASCII or UTF-8 already

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

commit 1171aa3961b5eb46a5d2ee64751c02a393a8a6f5
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Fri Jun 17 10:42:39 2016 +0200

correct short_ref document about include file

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

commit 2c6b2aac2ef04431ece2c6ff80e574aa2e58814b
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Fri Jun 17 10:40:34 2016 +0200

correct error message to new semantics

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

commit 506cde73451740bf0a15eff7d4abb158ee719ab0
Author: mflehmig <martin.flehmig@…>
Date: Fri Jun 17 10:14:26 2016 +0200

Fixed include of ADOL-C header.

ADOL-C header was included in old fashion (without adolc directory) for this example.

commit 2a023d3281d3d6d9824bad724a5768e3ee2fff94
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Thu Jun 16 13:50:39 2016 +0200

Try to use boost::pool for allocating advals in traceless vector mode

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

commit 80f1e2019ac1faab96fe06f3e9da47efcc1bcd23
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Mon May 23 15:13:22 2016 +0200

correct a reference in doc and rebuild

commit d7ab5283afe58bacb2e8739d72ede4e17f4c8081
Author: Mladen Banovic <mladenbanovic2705@…>
Date: Fri May 20 16:42:13 2016 +0200

Update section 7 of adolc-manual related to the Traceless forward differentiation.

commit bedb8e36f959c5272e4610fe504acc83208e5e9d
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Tue May 17 16:09:36 2016 +0200

macro name correction

commit 92ff596a0331776901df7f172ca347572e3daafd
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Tue May 17 15:56:17 2016 +0200

Add a warning about using static build of ADOL-C

static build of ADOL-C does not call constructors
for internal global objects, thereby causing
segmentation faults.

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

  • Property svn:keywords set to Author Date Id Revision
File size: 5.1 KB
Line 
1/*----------------------------------------------------------------------------
2 ADOL-C -- Automatic Differentiation by Overloading in C++
3 File:     liborpar.cpp
4 Revision: $Id: liborpar.cpp 708 2016-07-12 08:18:44Z kulshres $
5 Contents: example for differentiation of OpemMP parallel programs
6           parallel version
7
8 Copyright (c) Andrea Walther
9 
10 This file is part of ADOL-C. This software is provided as open source.
11 Any use, reproduction, or distribution of the software constitutes
12 recipient's acceptance of the terms of the accompanying license file.
13 
14---------------------------------------------------------------------------*/
15
16/* Program to compute deltas and vegas of swaption portfolio
17   from forward and reverse mode pathwise sensitivities
18   in parallel written by Andrea Walther in 2008-11 based on
19   code written by written by Mike Giles in 2005-7 which is
20   again based on code written by Zhao and Glasserman at
21   Columbia University  */
22
23using namespace std;
24
25#include <stdlib.h>
26#include <stdio.h>
27#include <ctime>
28#include <cmath>
29
30#include "adolc/adolc.h"
31
32#ifdef _OPENMP
33#include <omp.h>
34#include <adolc/adolc_openmp.h>
35#endif
36
37/* calculate path values */
38
39template <typename ADdouble>
40void path_calc(const int N, const int Nmat, const double delta,
41               ADdouble L[], const double lambda[], ADdouble z[])
42{
43  int      i, n;
44  double   lam, con1;
45  ADdouble v, vrat;
46  ADdouble sqez;
47 
48  for(n=0; n<Nmat; n++) {
49    sqez = sqrt(delta)*z[n];
50   
51    v = 0.0;
52    for (i=n+1; i<N; i++) {
53      lam  = lambda[i-n-1];
54      con1 = delta*lam;
55      v   += (con1*L[i])/(1.0+delta*L[i]);
56      vrat = exp(con1*v + lam*(sqez-0.5*con1));
57      L[i] = L[i]*vrat;
58    }
59  }
60}
61
62/* calculate the portfolio value v */
63
64template <typename ADdouble>
65void portfolio(const int N, const int Nmat, const double delta,
66               const int Nopt, const int maturities[],
67               const double swaprates[],
68               const ADdouble L[], ADdouble& v )
69{
70  int      i, m, n;
71  ADdouble b, s, swapval, *B, *S;
72  B = new ADdouble[N];
73  S = new ADdouble[N];
74 
75  b = 1.0;
76  s = 0.0;
77 
78  for (n=Nmat; n<N; n++) {
79    b    = b/(1.0+delta*L[n]);
80    s    = s + delta*b;
81    B[n] = b;
82    S[n] = s;
83  }
84 
85  v = 0;
86 
87  for (i=0; i<Nopt; i++){
88    m = maturities[i] + Nmat-1;
89    swapval = B[m] + swaprates[i]*S[m] - 1.0;
90    condassign(v,-swapval,v-100.0*swapval);
91  }
92 
93  // apply discount //
94 
95  for (n=0; n<Nmat; n++)
96    v = v/(1.0+delta*L[n]);
97  delete[](B); 
98  delete[](S); 
99}
100
101/* -------------------------------------------------------- */
102
103int main(){
104
105  // constants for all threads
106
107  // LIBOR interval  //
108  double delta = 0.25; 
109  // data for swaption portfolio //
110  int    Nopt = 15;
111  int    maturities[] = {4,4,4,8,8,8,20,20,20,28,28,28,40,40,40};
112  double swaprates[]  = {.045,.05,.055,.045,.05,.055,.045,.05,
113                         .055,.045,.05,.055,.045,.05,.055 };
114
115  int       i, j, N, Nmat, npath, nthreads, slot;
116  double    vtot, *v, *lambda, **z, **grad, *gradtot, **xp;
117
118  Nmat = 40;
119  N = Nmat+40;
120  npath = 30;
121  nthreads = 5;
122  slot = npath/nthreads;
123
124  lambda   = new double[N];
125  v        = new double[npath];
126  gradtot  = new double[N];
127  z        = new double*[npath];
128  grad     = new double*[npath];
129  xp       = new double*[npath];
130  for (i=0;i<npath;i++) 
131    {
132      z[i] = new double[N];
133      grad[i] = new double[N+Nmat];
134      xp[i] = new double[N+Nmat];
135   }
136
137  for (i=0;i<N;i++) 
138    {
139      gradtot[i] = 0.0;
140      lambda[i] = 0.2;
141    }
142
143  for (j=0; j<npath; j++)
144    {
145      v[j] = 0;
146      for (i=0; i<N; i++) 
147          xp[j][i]=  0.05;
148      for (i=0; i<Nmat; i++) 
149        {
150          z[j][i] = 0.3;
151          xp[j][N+i]=  0.3;
152        }
153    }
154
155  omp_set_num_threads(nthreads);
156
157  //----------------------------------------------------------//
158  //                                                          //
159  // do a full path + portfolio sensitivity check             //
160  //                                                          //
161  // A real application would generate a different random     //
162  // vector z for each path but here we set one and reuse it  //
163  //                                                          //
164  //----------------------------------------------------------//
165
166#pragma omp parallel firstprivate(ADOLC_OpenMP_Handler)
167  {   
168    // different paths for each thread
169    int index = omp_get_thread_num();
170    int l,k;
171    int rv = 0;
172   
173    adouble *La, va, *za;
174
175    La = new adouble[N];
176    za  = new adouble[Nmat];
177
178    int init = index*slot;
179
180    trace_on(init);
181      for(k=0;k<N;k++) 
182        La[k] <<= 0.050000;
183      for(j=0;j<Nmat;j++) 
184        za[j] <<= z[init][j];
185   
186      path_calc(N,Nmat,delta,La,lambda,za);
187      portfolio(N,Nmat,delta,Nopt,maturities,swaprates,La,va);
188       
189      va >>= v[init];
190    trace_off(1);
191   
192    for(l=init;l<init+slot;l++)
193      {
194        rv = gradient(init,N+Nmat,xp[l],grad[l]);
195      }
196
197#pragma omp barrier
198  if (index==0)
199    {
200      vtot = 0;
201      for (l=0; l<npath; l++)
202        {
203          vtot += v[l];
204          for(k=0;k<N;k++)
205            gradtot[k] += grad[l][k];
206        }
207      vtot = vtot/npath;
208      for(k=0;k<N;k++)
209        gradtot[k] /= npath;
210
211      printf("Gradient: \n");
212      for(i=0;i<N;i++)
213        printf(" %f \n",gradtot[i]);
214    }
215  } 
216
217
218  return 0;
219}
220
Note: See TracBrowser for help on using the repository browser.