source: trunk/ADOL-C/examples/additional_examples/timing/vfunc_gear.cpp @ 171

Last change on this file since 171 was 171, checked in by kulshres, 9 years ago

Squashed merge branch 'master' of 'gitclone' into svn

  • 'master' of 'gitclone': (84 commits) adjust example makefiles and include paths get rid of the symlink in the src subdirectory

details of the commits:
commit c9e4bc332d2363f737fc2e8a8fcfc2e43ddb9d15
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Mon Oct 4 15:43:47 2010 +0200

adjust example makefiles and include paths

include paths in example sources were wrong for some time now
simplify makefile rules too, there is really no need for checking SPARSE
adjust include paths in makefiles.

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

commit e6e1963e41e097fd5b4a79cd1611c12f6868dc94
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Mon Oct 4 15:41:25 2010 +0200

get rid of the symlink in the src subdirectory

windows doesn't like symlinks and make infinite depth directories
we now create a symlink for build in the directory parallel to src
adjust all makefiles.am accordingly for build

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

  • Property svn:keywords set to Author Date Id Revision
File size: 11.3 KB
Line 
1/*----------------------------------------------------------------------------
2 ADOL-C -- Automatic Differentiation by Overloading in C++
3 File:     vfunc_gear.cpp
4 Revision: $Id: vfunc_gear.cpp 171 2010-10-04 13:57:19Z kulshres $
5 Contents: Example of function module containing the machine tool example
6            of gearing
7 
8   Each << function module >> contains:
9         
10     (1) const char* const controlFileName
11     (2) int indepDim;
12     (3) int depDim;
13     (4) void initProblemParameters( void )
14     (5) void initIndependents( double* indeps )
15     (6) void originalVectorFunction( double* indeps, double* deps )
16     (7) void tapingVectorFunction( int tag, double* indeps, double* deps )   
17 
18 Copyright (c) Andrea Walther, Andreas Griewank, Andreas Kowarz,
19               Hristo Mitev, Sebastian Schlenkrich, Jean Utke, Olaf Vogel
20
21 This file is part of ADOL-C. This software is provided as open source.
22 Any use, reproduction, or distribution of the software constitutes
23 recipient's acceptance of the terms of the accompanying license file.
24 
25---------------------------------------------------------------------------*/
26#define _VFUNC_GEAR_C_
27
28/****************************************************************************/
29/*                                                                 INCLUDES */
30#include <adolc/adolc.h>
31
32#include <cmath>
33#include <time.h>
34#include <cstdlib>
35
36
37/****************************************************************************/
38/*                                                         GLOBAL VARIABLES */
39
40/*--------------------------------------------------------------------------*/
41/*                                                        Control file name */
42const char* controlFileName = "gearexam.ctrl";
43
44/*--------------------------------------------------------------------------*/
45/*                                                               Dimensions */
46int indepDim;
47int depDim;
48
49/*--------------------------------------------------------------------------*/
50/*                                       Other problem dependent parameters */
51//static unsigned short int dx[3]; /* variable needed by erand48(.) */
52
53
54/****************************************************************************/
55/*                                                  INIT PROBLEM PARAMETERS */
56void initProblemParameters( void ) {
57    fprintf(stdout,"GEAREXAM (ADOL-C Example)\n\n");
58
59    /* number of indeps & deps */
60    indepDim = 3;
61    depDim   = 3;
62
63    /* Init erand48(); */
64    struct tm s;
65    time_t t;
66    time(&t);
67    s=*localtime(&t);
68    srand(s.tm_sec*s.tm_min);
69    /*  dx[0]=rand();
70      dx[1]=rand();
71      dx[2]=rand();*/
72}
73
74
75/****************************************************************************/
76/*                                                        INITIALIZE INDEPs */
77void initIndependents( double* indeps ) {
78    for (int i=0; i<indepDim; i++)
79        indeps[i] = (double)rand();
80}
81
82
83/****************************************************************************/
84/*                                                 ORIGINAL SCALAR FUNCTION */
85
86/*--------------------------------------------------------------------------*/
87/*                                                       The model function */
88#define Pi 3.141592654
89
90/*--------------------------------------------------------------------------*/
91// important machine tool parameters
92
93int konvex=     1;          // konvexe oder konkave Flanke
94
95int zz=         13;         // Zaehnezahl
96double delta=   16.46;      // Kegelwinkel
97double st=      0.00001;    // Verschiebung Teilkegelspitze
98
99double xmw=     79.33279;   // MK-x
100double ymw=     -74.05;     // MK-y
101double zmw=     159.49741;  // MK-z
102double m=       0.0;        // Erzeugungs-Achsversatz
103double zwr=     0.0;        // Reitstockeinstellung
104double theta_w= -58.2253;   // Waelztrommelwinkel
105
106double xmk=     17.49874;   // Messerversatz
107int zm=         5;          // MK-Gangzahl
108double ymk=     80.0;       // MK-Versatz
109
110double rho=     101.64155;  // Spitzenradius, Flugkreisr.
111double r=       2.1;        // Kopfradius
112double rs=      2594.425;   // Sphaerikradius
113double ys=      876.147;    // Sphaerik-Mitte-Y
114double zs=      -2442.015;  // Sphaerik-Mitte-Z
115
116
117/*--------------------------------------------------------------------------*/
118// elementary rotations
119
120void D1 ( adouble * vec, double & alpha ) {
121    double locCos=cos(alpha);
122    double locSin=sin(alpha);
123    adouble tmpVec2=locSin*vec[1] + locCos*vec[2];
124    vec[1]=locCos*vec[1] - locSin*vec[2];
125    vec[2]=tmpVec2;
126}
127
128void D2 ( adouble * vec, double & alpha ) {
129    double locCos=cos(alpha);
130    double locSin=sin(alpha);
131    adouble tmpVec2=-locSin*vec[0] + locCos*vec[2];
132    vec[0]=locCos*vec[0] + locSin*vec[2];
133    vec[2]=tmpVec2;
134}
135
136void D3 ( adouble * vec, double & alpha ) {
137    double locCos=cos(alpha);
138    double locSin=sin(alpha);
139    adouble tmpVec1=locSin*vec[0] + locCos*vec[1];
140    vec[0]=locCos*vec[0] - locSin*vec[1];
141    vec[1]=tmpVec1;
142}
143
144void D1 ( adouble * vec, adouble & alpha ) {
145    adouble locCos=cos(alpha);
146    adouble locSin=sin(alpha);
147    adouble tmpVec2=locSin*vec[1] + locCos*vec[2];
148    vec[1]=locCos*vec[1] - locSin*vec[2];
149    vec[2]=tmpVec2;
150}
151
152void D2 ( adouble * vec, adouble & alpha ) {
153    adouble locCos=cos(alpha);
154    adouble locSin=sin(alpha);
155    adouble tmpVec2=-locSin*vec[0] + locCos*vec[2];
156    vec[0]=locCos*vec[0] + locSin*vec[2];
157    vec[2]=tmpVec2;
158}
159
160void D2 ( adouble *  depVec, adouble * indepVec,  adouble & alpha ) {
161    if ( indepVec == depVec ) {
162        D2(depVec,alpha);
163        return;
164    }
165    adouble locCos=cos(alpha);
166    adouble locSin=sin(alpha);
167    depVec[0]=locCos*indepVec[0] + locSin*indepVec[2];
168    depVec[1]=indepVec[1];
169    depVec[2]=-locSin*indepVec[0] + locCos*indepVec[2];
170}
171
172void D3 ( adouble * vec, adouble & alpha ) {
173    adouble locCos=cos(alpha);
174    adouble locSin=sin(alpha);
175    adouble tmpVec1=locSin*vec[0] + locCos*vec[1];
176    vec[0]=locCos*vec[0] - locSin*vec[1];
177    vec[1]=tmpVec1;
178}
179
180void D1 ( double * vec, double & alpha ) {
181    double locCos=cos(alpha);
182    double locSin=sin(alpha);
183    double tmpVec2=locSin*vec[1] + locCos*vec[2];
184    vec[1]=locCos*vec[1] - locSin*vec[2];
185    vec[2]=tmpVec2;
186}
187
188void D2 ( double * vec, double & alpha ) {
189    double locCos=cos(alpha);
190    double locSin=sin(alpha);
191    double tmpVec2=-locSin*vec[0] + locCos*vec[2];
192    vec[0]=locCos*vec[0] + locSin*vec[2];
193    vec[2]=tmpVec2;
194}
195
196void D2 ( double * depVec, double * indepVec,  double & alpha ) {
197    if ( indepVec == depVec ) {
198        D2(depVec,alpha);
199        return;
200    }
201    double locCos=cos(alpha);
202    double locSin=sin(alpha);
203    depVec[0]=locCos*indepVec[0] + locSin*indepVec[2];
204    depVec[1]=indepVec[1];
205    depVec[2]=-locSin*indepVec[0] + locCos*indepVec[2];
206}
207
208void D3 ( double * vec, double & alpha ) {
209    double locCos=cos(alpha);
210    double locSin=sin(alpha);
211    double tmpVec1=locSin*vec[0] + locCos*vec[1];
212    vec[0]=locCos*vec[0] - locSin*vec[1];
213    vec[1]=tmpVec1;
214}
215
216
217/*--------------------------------------------------------------------------*/
218// parametrized cutting edge
219
220void def_messer(adouble *z, adouble *messer) {
221
222    double u0, uOri, phi0;
223    adouble h;
224
225    phi0= asin((r+ys)/(r+rs));
226    if (konvex==1) {
227        u0=rs*phi0;
228        uOri=1.0;
229    } else {
230        u0=rs*(phi0-Pi);
231        uOri=-1.0;
232    };
233
234    h= (z[0]-u0)/(uOri*rs);
235    messer[0]=zs+rs*cos(h);
236    messer[1]=0.0;
237    messer[2]=-ys-rs*sin(h);
238}
239
240void def_messer(double *z, double *messer) {
241
242    double u0, uOri, phi0;
243    double h;
244
245    phi0= asin((r+ys)/(r+rs));
246    if (konvex==1) {
247        u0=rs*phi0;
248        uOri=1.0;
249    } else {
250        u0=rs*(phi0-Pi);
251        uOri=-1.0;
252    };
253
254    h= (z[0]-u0)/(uOri*rs);
255    messer[0]=zs+rs*cos(h);
256    messer[1]=0.0;
257    messer[2]=-ys-rs*sin(h);
258}
259
260/*--------------------------------------------------------------------------*/
261// the main function
262
263void gearFunction(double* pz, double* pf) {
264    double ah;
265    double* messer;
266    messer= new double[3];
267
268    def_messer(pz, messer);
269
270    // Position der Schneide am Messerkopf
271    messer[0]+=rho;
272    messer[1]-=xmk;
273
274    // Messerkopfrotation mit Parameter v
275    D3(messer,pz[1]);
276
277    // Lage des Messerkopfs auf der Wiege
278    messer[2]-=ymk;
279
280    // Eindrehen in Orientierung der Wiege
281    ah=messer[0];
282    messer[0]=messer[1];
283    messer[1]=ah;
284    messer[2]=-messer[2];
285
286    // Verschiebung
287    messer[0]-=xmw;
288    messer[1]+=zmw;
289    messer[2]+=ymw;
290
291    // Wiegenwinkel thetaW, entspricht dem wert t=0
292    // + Wiegenbewegung mit Parameter t
293    ah = theta_w+pz[2];
294    D3(messer,ah);
295
296    // Achsversatz
297    messer[0]+=m;
298
299    // Eindrehen in Orientierung des Werkrades
300    messer[0]=-messer[0];
301    ah=messer[1];
302    messer[1]=-messer[2];
303    messer[2]=-ah;
304
305    // Teilkegeloeffnungswinkel delta, y-Achsen entgegengesetzt
306    D1(messer,delta);
307
308    // neue Verschiebung der Werkradachse
309    messer[2]+=zwr+st;
310
311    // gekoppelte Werkraddrehung in Abhaengigkeit von t und v
312    ah = (1/sin(delta))*pz[2] + (double(zm)/zz)*pz[1];
313    D2(pf,messer,ah);
314
315    delete[] messer;
316}
317
318/*--------------------------------------------------------------------------*/
319/*                                                   The interface function */
320void originalVectorFunction( double* indeps, double* deps ) {
321    gearFunction(indeps,deps);
322}
323
324
325/****************************************************************************/
326/*                                                   TAPING SCALAR FUNCTION */
327
328/*--------------------------------------------------------------------------*/
329/*                                                       The model function */
330void activeGearFunction(adouble* z, adouble* f) {
331    adouble ah;
332
333    adouble* messer = new adouble[3];
334    def_messer(z, messer);
335
336    // Position der Schneide am Messerkopf
337    messer[0]+=rho;
338    messer[1]-=xmk;
339
340    // Messerkopfrotation mit Parameter v
341    D3(messer,z[1]);
342
343    // Lage des Messerkopfs auf der Wiege
344    messer[2]-=ymk;
345
346    // Eindrehen in Orientierung der Wiege
347    ah=messer[0];
348    messer[0]=messer[1];
349    messer[1]=ah;
350    messer[2]=-messer[2];
351
352    // Verschiebung
353    messer[0]-=xmw;
354    messer[1]+=zmw;
355    messer[2]+=ymw;
356
357    // Wiegenwinkel thetaW, entspricht dem wert t=0
358    // + Wiegenbewegung mit Parameter t
359    ah = theta_w+z[2];
360    D3(messer,ah);
361
362    // Achsversatz
363    messer[0]+=m;
364
365    // Eindrehen in Orientierung des Werkrades
366    messer[0]=-messer[0];
367    ah=messer[1];
368    messer[1]=-messer[2];
369    messer[2]=-ah;
370
371    // Teilkegeloeffnungswinkel delta, y-Achsen entgegengesetzt
372    D1(messer,delta);
373
374    // neue Verschiebung der Werkradachse
375    messer[2]+=zwr+st;
376
377    // gekoppelte Werkraddrehung in Abhaengigkeit von t und v
378    ah = (1/sin(delta))*z[2] + (double(zm)/zz)*z[1];
379    D2(f,messer,ah);
380
381    delete[] messer;
382}
383
384
385/*--------------------------------------------------------------------------*/
386/*                                                   The interface function */
387void tapingVectorFunction( int tag, double* indeps, double* deps ) {
388    int i;
389    trace_on(tag);
390    adouble* activeIndeps = new adouble[indepDim];
391    adouble* activeDeps   = new adouble[depDim];
392    adouble* aIP = activeIndeps;
393    double*  iP  = indeps;
394    for (i=0; i<indepDim; i++)
395        *aIP++ <<= *iP++;
396    activeGearFunction(activeIndeps,activeDeps);
397    aIP = activeDeps;
398    iP  = deps;
399    for (i=0; i<depDim; i++)
400        *aIP++ >>= *iP++;
401    trace_off();
402}
403
404#undef _VFUNC_GEAR_C_
405
406
407
408
409
Note: See TracBrowser for help on using the repository browser.