1 """Variable is a one-dimensional discrete and continuous real variable class.
2
3 Robert Clewley, July 2005
4 """
5
6
7
8
9 from utils import *
10 from common import *
11 from errors import *
12 from Points import *
13 from Interval import *
14 from FuncSpec import ImpFuncSpec
15
16 from numpy import Inf, NaN, isfinite, sometrue, alltrue, any, all, \
17 array, float64, int32, ndarray, asarray
18
19 import copy
20 import types, math, random
21
22 __all__ = ['Variable', 'HybridVariable',
23 'OutputFn', 'isinputcts', 'isinputdiscrete',
24 'isoutputcts', 'isoutputdiscrete',
25 'iscontinuous', 'isdiscrete',
26 'numeric_to_vars', 'pointset_to_vars']
27
28
29
30
33 if self.warnings != []:
34 output = "Warnings:"
35 for (i,d) in self.warnings:
36 if d is None:
37 output += "Independent variable value %s was out of "% i + \
38 "bounds"
39 else:
40 output += "Dependent variable value %s was out of " % s + \
41 "bounds at independent variable value %s" % i
42 else:
43 output = ''
44 return output
45
46
48 """Utility to convert Pointset to a dictionary of Variables.
49 If discrete option set to False (default is True) then the
50 Variables will be linearly interpolated within their domain.
51
52 Any labels in the pointset will be preserved in the Variables
53 in case of their re-extraction using the getDataPoints method.
54 """
55 coordnames = pts.coordnames
56 vals = pts.coordarray
57 all_types_float = pts.coordtype == float
58 if isparameterized(pts):
59 indepvar = pts.indepvararray
60 indepvarname = pts.indepvarname
61 if discrete:
62 indepvartype = int
63 else:
64 indepvartype = float
65 indepdomain = Interval(pts.indepvarname, indepvartype,
66 extent(pts.indepvararray),
67 abseps=pts._abseps)
68 else:
69 indepvar = None
70 indepvarname = None
71 indepdomain = None
72 return numeric_to_vars(vals, coordnames, indepvar, indepvarname,
73 indepdomain, all_types_float, discrete,
74 pts._abseps, pts.labels)
75
76
77 -def numeric_to_vars(vals, coordnames, indepvar=None, indepvarname='t',
78 indepdomain=None, all_types_float=True, discrete=True,
79 abseps=None, labels=None):
80 """Utility to convert numeric types to a dictionary of Variables.
81 If discrete option set to True (default is False) then the
82 Variables will be linearly interpolated within their domain.
83 """
84 if isinstance(coordnames, str):
85 coordnames = [coordnames]
86 if isinstance(vals, _num_types):
87 vals = [[vals]]
88 vars = {}
89 if indepvar is None:
90 for i, c in enumerate(coordnames):
91 if all_types_float:
92 vartype = float
93 else:
94 vartype = array(vals[i]).dtype.type
95 if discrete:
96 vars[c] = Variable(outputdata=Pointset({'coordnames': c,
97 'coordarray': vals[i],
98 'coordtype': vartype}),
99 name=c, abseps=abseps,
100 labels=labels)
101 else:
102 raise AssertionError("Cannot use continuously defined "
103 "option without an independent variable")
104 return vars
105 else:
106 if isinstance(indepvar, _num_types):
107 indepvartype = type(indepvar)
108 indepvar = [indepvar]
109 else:
110 indepvartype = asarray(indepvar).dtype.type
111 if indepdomain is None:
112 indepdomain = indepvarname
113 else:
114 if isinstance(indepdomain, Interval):
115 assert indepvarname == indepdomain.name, "Indep varname mismatch"
116 else:
117 if discrete:
118 var_type = int
119 else:
120 var_type = float
121 indepdomain = Interval(indepvarname, var_type, indepdomain)
122 for i, c in enumerate(coordnames):
123 if all_types_float:
124 vartype = float
125 else:
126 vartype = array(vals[i]).dtype.type
127 if discrete:
128 vars[c] = Variable(outputdata=Pointset({'coordnames': c,
129 'coordarray': vals[i],
130 'coordtype': vartype,
131 'indepvarname': indepvarname,
132 'indepvararray': indepvar,
133 'indepvartype': indepvartype}),
134 indepdomain=indepdomain, name=c,
135 abseps=abseps, labels=labels)
136 else:
137 dom_int = Interval(c, vartype, extent(vals[i]),
138 abseps=abseps)
139 vars[c] = Variable(outputdata=interp1d(indepvar, vals[i]),
140 indepdomain=indepdomain,
141 depdomain=dom_int, name=c,
142 abseps=abseps, labels=labels)
143 return vars
144
145
147 """One-dimensional discrete and continuous real variable class.
148 """
149
150 - def __init__(self, outputdata=None, indepdomain=None, depdomain=None,
151 name='noname', abseps=None, labels=None):
152
153
154 self._funcreg = {}
155 if isinstance(name, str):
156
157 self.name = name
158 else:
159 raise TypeError("name argument must be a string")
160
161 if outputdata is None or isinstance(outputdata, (Pointset, interp1d)):
162 if indepdomain is None:
163 indepdomain = 't'
164 if depdomain is None:
165 depdomain = 'x'
166
167
168 self._vectorizable = True
169 self.defined = False
170 self.indepdomain = None
171 self.indepvartype = None
172 self.indepvarname = None
173 self.depdomain = None
174 self.coordtype = None
175 self.coordname = None
176 self._refvars = None
177
178 self.trajirange = None
179 self.trajdrange = None
180
181 self.setIndepdomain(indepdomain, abseps)
182
183 self._internal_t_offset = 0
184
185 self.setOutput(outputdata, abseps)
186
187 self.setDepdomain(depdomain, abseps)
188 assert self.coordname != self.indepvarname, ("Independent variable "
189 "name and coordinate name must be different")
190 self.diagnostics = VarDiagnostics()
191
192
193
194 self.labels = labels
195
196
199
202
203
204
206 return self.globalt0 + t
207
208
210 return self.initialconditions[varname]
211
212
214 if x>0:
215 return 1
216 else:
217 return 0
218
219
221 if c:
222 return e1
223 else:
224 return e2
225
226
228 return self._var_namemap[varname]
229
230
232 """Add dynamically-created methods to Veriable object"""
233
234
235 for auxfnname in funcspec.auxfns:
236 fninfo = funcspec.auxfns[auxfnname]
237 if not hasattr(Variable, fninfo[1]):
238
239
240 try:
241 exec fninfo[0] in globals()
242 except:
243 print 'Error in supplied auxiliary function code'
244 raise
245 self._funcreg[fninfo[1]] = ('Variable', fninfo[0])
246 setattr(Variable, fninfo[1], eval(fninfo[1]))
247
248 fninfo_spec = funcspec.spec
249 if not hasattr(Variable, fninfo_spec[1]):
250 try:
251 exec fninfo_spec[0] in globals()
252 except:
253 print 'Error in supplied functional specification code'
254 raise
255 self._funcreg[fninfo_spec[1]] = ('Variable', fninfo_spec[0])
256 setattr(Variable, fninfo_spec[1], eval(fninfo_spec[1]))
257
258 if funcspec.auxspec:
259 fninfo_auxspec = funcspec.auxspec
260 if not hasattr(Variable, fninfo_auxspec[1]):
261 try:
262 exec fninfo_auxspec[0] in globals()
263 except:
264 print 'Error in supplied auxiliary variable code'
265 raise
266 self._funcreg[fninfo_auxspec[1]] = ('Variable', fninfo_auxspec[0])
267 setattr(Variable, fninfo_auxspec[1], eval(fninfo_auxspec[1]))
268
269 if isinstance(funcspec, ImpFuncSpec):
270 impfn_name = funcspec.algparams['impfn_name']
271 if funcspec.algparams['jac']:
272 jac_str = "fprime=funcspec.algparams['jac'],"
273 else:
274 jac_str = ""
275
276
277
278
279
280
281 if len(funcspec.vars) == 1:
282
283
284
285
286 specfn_str = "lambda a1, a2, a3: " \
287 + fninfo_spec[1] \
288 + "(None, a2, [a1], a3)[0]"
289 else:
290
291 specfn_str = "lambda a1, a2, a3: " \
292 + fninfo_spec[1] \
293 + "(None, a2, a1, a3)"
294 this_scope = globals()
295 this_scope.update({'funcspec': locals()['funcspec'],
296 'fninfo_spec': locals()['fninfo_spec']})
297 impfn_str = impfn_name + \
298 " = makeImplicitFunc(" + specfn_str + "," \
299 + jac_str + """x0=funcspec.algparams['x0'],
300 extrafargs=(funcspec.algparams['pars'],),
301 xtolval=funcspec.algparams['atol'],
302 maxnumiter=funcspec.algparams['maxnumiter'],
303 solmethod=funcspec.algparams['solvemethod'],
304 standalone=False)"""
305 try:
306 exec impfn_str in this_scope
307 except:
308 print 'Error in supplied implicit function code'
309 raise
310
311
312 self._funcreg['_impfn'] = (impfn_name, impfn_str)
313
314 setattr(Variable, impfn_name, eval(impfn_name))
315
316 del this_scope['funcspec']
317 del this_scope['fninfo_spec']
318
320 """Reveal underlying mesh and values at mesh points, provided
321 Variable is based on a mesh (otherwise None is returned).
322 The returned Pointset will be time-shifted according to the
323 Variable's current _internal_t_offset attribute.
324
325 Any pointset labels present when the variable was created will
326 be restored.
327 """
328 if isinstance(self.output, VarCaller):
329 return Pointset(indepvarname=self.indepvarname,
330 indepvararray=self.output.pts.indepvararray + self._internal_t_offset,
331 coordnames=[self.coordname],
332 coordarray=self.output.pts.coordarray[0],
333 labels=self.labels)
334 elif hasattr(self.output, 'datapoints'):
335 datapoints = self.output.datapoints
336 return Pointset(indepvarname=self.indepvarname,
337 indepvararray=datapoints[0] + self._internal_t_offset,
338 coordnames=[self.coordname],
339 coordarray=datapoints[1],
340 labels=self.labels)
341 else:
342 return None
343
345 """Reveal underlying mesh as arrays, rather than Pointset
346 as returned by getDataPoints method. If no underlying mesh is
347 present, None is returned."""
348 try:
349
350 mesh = self.output.datapoints
351 except AttributeError:
352 try:
353
354 pts = self.output.pts
355 mesh = array([pts.indepvararray, pts.coordarray[0]])
356 except AttributeError:
357 mesh = None
358 return mesh
359
361 mesh = self.underlyingMesh()
362 if mesh is None:
363 raise RuntimeError("Cannot truncate a Variable without an underlying mesh by index")
364 try:
365 new_t_end = mesh[0][idx]
366 except IndexError:
367 raise ValueError("Truncation index %d out of range"%idx)
368 except TypeError:
369 raise TypeError("Index must be an integer")
370 if isinstance(self.indepdomain, Interval):
371 self.indepdomain.set([self.indepdomain[0], new_t_end])
372 else:
373
374 self.indepdomain = self.indepdomain[0:idx]
375
376
377 if isinstance(self.depdomain, ndarray):
378 self.depdomain = self.depdomain[0:idx]
379
380 self._setRanges(self.depdomain._abseps)
381
383
384
385 try:
386 output = self.output
387 except AttributeError:
388
389 return
390 if isinstance(output, VarCaller):
391 self.trajirange = Interval('traj_indep_bd',
392 self.indepvartype,
393 extent(output.pts.indepvararray),
394 abseps=abseps)
395 self.trajdrange = Interval('traj_dep_bd',
396 self.coordtype,
397 extent(output.pts.coordarray[0]),
398 abseps=abseps)
399 elif type(output) in [types.InstanceType, types.TypeType] \
400 or isinstance(output, (OutputFn, interpclass)):
401 if hasattr(output, 'types'):
402 deptype = output.types[0]
403 indeptype = output.types[1]
404 else:
405
406 deptype = indeptype = float
407 if isinstance(output.datapoints[0], Interval):
408 assert compareNumTypes(output.types[0], \
409 output.datapoints[0].type), \
410 "Inconsistent type with Interval bounds"
411 self.trajirange = output.datapoints[0]
412 else:
413 self.trajirange = Interval('traj_indep_bd', indeptype,
414 extent(output.datapoints[0]),
415 abseps=abseps)
416 if isinstance(output.datapoints[1], Interval):
417 assert compareNumTypes(output.types[1], \
418 output.datapoints[1].type), \
419 "Inconsistent type with Interval bounds"
420 self.trajdrange = output.datapoints[1]
421 else:
422 self.trajdrange = Interval('traj_dep_bd', deptype,
423 extent(output.datapoints[1]),
424 abseps=abseps)
425
426 - def setOutput(self, outputdata, funcspec=None, globalt0=0,
427 var_namemap=None, ics=None, refvars=None, abseps=None):
428 """Dynamically create 'output' method of Variable"""
429
430 self.globalt0 = globalt0
431 if type(outputdata) in [types.FunctionType,
432 types.BuiltinFunctionType,
433 types.MethodType]:
434
435 self.output = outputdata
436 assert ics is None, "Invalid option for this type of output"
437 if outputdata != noneFn:
438 self.defined = True
439 elif isinstance(outputdata, tuple):
440
441
442 assert len(outputdata) == 2, "Incorrect size of outputdata tuple"
443 if funcspec is not None:
444 self.addMethods(funcspec)
445 self._var_namemap = var_namemap
446 self._funcreg['funcspec'] = (None, funcspec)
447 else:
448 raise ValueError('funcspec missing in setOutput')
449
450 try:
451 exec outputdata[1] in globals()
452 except:
453 print 'Internal Error in _mapspecfn code'
454 raise
455 has_op = hasattr(self, 'output')
456
457
458 if not has_op or (has_op and self.output is noneFn):
459 def wrap_output(arg):
460 return eval(outputdata[0])(self, arg)
461 setattr(self, 'output', wrap_output)
462 self._funcreg['outputdata'] = (None, outputdata)
463 t0 = self.indepdomain[0]
464 if ics is None and not isinstance(funcspec, ImpFuncSpec):
465 try:
466 self.initialconditions = {self.coordname: self.output(t0)}
467 except ValueError:
468 self.initialconditions = {self.coordname: NaN}
469 except TypeError:
470 print "Debugging info: self.output = ", self.output
471 raise
472 else:
473 self.initialconditions = ics
474 self._vectorizable = False
475 self._refvars = refvars
476 self.defined = True
477 elif type(outputdata) in [types.InstanceType, types.TypeType] \
478 or isinstance(outputdata, (OutputFn, interpclass)):
479
480
481
482
483 assert ics is None, "Invalid option for this type of output"
484 assert '__call__' in dir(outputdata), "Must provide callable object"
485 self.output = outputdata
486 if hasattr(outputdata, 'datapoints'):
487 self._setRanges(abseps)
488 self.defined = True
489 elif isinstance(outputdata, Pointset):
490
491 assert ics is None, "Invalid option for this type of output"
492 assert isparameterized(outputdata), ("Must only pass parameterized"
493 " pointsets")
494 if outputdata.dimension == 1:
495 self.coordname = copy.copy(outputdata.coordnames[0])
496 self.indepvarname = outputdata.indepvarname
497 self.output = VarCaller(outputdata)
498 self.coordtype = outputdata.coordtype
499 self.indepvartype = outputdata.indepvartype
500 if self.indepdomain is not None:
501 for v in outputdata[self.indepvarname]:
502 if not v in self.indepdomain:
503 raise ValueError("New Pointset data violates "
504 "independent variable domain already specified")
505 if self.depdomain is not None:
506 for v in outputdata[self.coordname]:
507 if not v in self.depdomain:
508 raise ValueError("New Pointset data violates "
509 "dependent variable domain already specified")
510 self._setRanges(abseps)
511 self.defined = True
512 else:
513 raise ValueError("Pointset data must be 1D to create a "
514 "Variable")
515 elif outputdata is None:
516
517 assert ics is None, "Invalid option when outputdata argument is None"
518 self.output = noneFn
519 self.defined = False
520 else:
521 raise TypeError("Invalid type for data argument: " \
522 +str(type(outputdata)))
523
524
525 - def setIndepdomain(self, indepdomain, abseps=None):
526 if isinstance(indepdomain, str):
527 self.indepvarname = indepdomain
528 if self.indepdomain is not None:
529
530
531 self.indepvarname = indepdomain
532 self.indepdomain.name = indepdomain
533 else:
534 self.indepdomain = Interval(self.indepvarname, float,
535 [-Inf, Inf], abseps=abseps)
536 self.indepvartype = float
537 else:
538 if isinstance(indepdomain, Interval):
539 if self.trajirange:
540 if indepdomain.contains(self.trajirange) is notcontained:
541 raise ValueError("Cannot set independent variable"
542 " domain inside current trajectory's"
543 " range")
544 self.indepdomain = indepdomain
545 self.indepvarname = indepdomain.name
546 self.indepvartype = _num_name2type[indepdomain.typestr]
547 elif isinstance(indepdomain, dict):
548
549 assert len(indepdomain) == 1, "Independent variable " \
550 "dictionary must have only 1 entry"
551 d = indepdomain.values()[0]
552 assert all(isfinite(d)), "Independent variable values must be" \
553 " finite"
554 if self.trajirange:
555 assert self.trajirange[0] in d
556 assert self.trajirange[1] in d
557 self.indepvarname = indepdomain.keys()[0]
558 if isinstance(d, (list, tuple)):
559 if self.coordtype is not None:
560 self.indepdomain = array(d, self.coordtype)
561 else:
562 self.indepdomain = array(d)
563 elif isinstance(d, ndarray):
564 da = array(d)
565 if self.indepvartype is not None and \
566 self.indepvartype != da.dtype.type:
567 raise TypeError("Mismatch between type of indepdomain "
568 "argument and Pointset data")
569 else:
570 self.indepdomain = da
571 else:
572 raise TypeError("Invalid type for independent "
573 "variable domain")
574
575
576
577 assert isincreasing(self.indepdomain), \
578 "Independent variable values must be increasing"
579 self.indepvartype = self.indepdomain.dtype.type
580 else:
581 print "Independent variable argument domain was:", indepdomain
582 raise TypeError("Invalid type for independent variable "
583 "domain")
584
585
586 - def setDepdomain(self, depdomain, abseps=None):
587 if isinstance(depdomain, str):
588 self.coordname = depdomain
589 if self.depdomain is None:
590 if self.coordtype is None:
591 self.depdomain = Interval(self.coordname, float,
592 [-Inf, Inf], abseps=abseps)
593 self.coordtype = float
594 else:
595 self.depdomain = Interval(self.coordname,
596 self.coordtype,
597 _num_maxmin[self.coordtype],
598 abseps=abseps)
599 else:
600
601
602 if isinstance(self.output, interpclass) and \
603 isinstance(self.depdomain, Interval):
604 self.depdomain.name = depdomain
605 else:
606 assert isinstance(self.output, Pointset)
607 self.diagnostics.warnings.append((self.depdomain.name,
608 "Dependent variable already named. "
609 "Ignoring user-supplied name."))
610 else:
611 if isinstance(depdomain, Interval):
612 if self.trajdrange:
613 if depdomain.contains(self.trajdrange) is notcontained:
614 raise ValueError("Cannot set dependent variable "
615 "domain inside current trajectory's "
616 "range")
617 self.depdomain = depdomain
618 self.coordname = depdomain.name
619 if self.coordtype is None:
620 self.coordtype = depdomain.type
621 elif self.coordtype == depdomain.type:
622 pass
623 else:
624 raise TypeError("Mismatch between type of depdomain "
625 "argument and Pointset coord data")
626 elif isinstance(depdomain, dict):
627 assert len(depdomain) == 1, \
628 "Depend variables dictionary must have only 1 entry"
629 d = depdomain.values()[0]
630 if self.trajdrange:
631 assert self.trajdrange[0] in d
632 assert self.trajdrange[1] in d
633
634 assert all(isfinite(d)), "Values must be finite"
635 self.coordname = depdomain.keys()[0]
636 if isinstance(d, (list, tuple)):
637 if self.coordtype is not None:
638 self.depdomain = array(d, self.coordtype)
639 else:
640 self.depdomain = array(d)
641 elif isinstance(d, ndarray):
642 da = array(d)
643 if self.coordtype is not None and \
644 self.coordtype != da.dtype.type:
645 raise TypeError("Mismatch between type of depdomain "
646 "argument and Pointset coord data")
647 else:
648 self.depdomain = da
649 else:
650 raise TypeError("Invalid type for dependent variable "
651 "domain")
652 self.coordtype = self.depdomain.dtype.type
653 else:
654 print "Dependent variable domain argument was:", depdomain
655 raise TypeError("Invalid type for dependent variable domain")
656 if isinstance(self.output, Pointset):
657 assert self.coordname == self.output.coordnames[0], \
658 "Mismatch between Pointset coord name and declared name"
659 assert self.indepvarname == self.output.indepvarname, \
660 ("Mismatch between Pointset independent variable name "
661 "and declared name")
662
663
664 - def __call__(self, indepvar, checklevel=0):
665
666
667
668
669 indepvar = asarray(indepvar) - self._internal_t_offset
670 if checklevel == 0:
671
672
673
674 try:
675 if not self._vectorizable and isinstance(indepvar, _seq_types):
676 return [self.output(ival) for ival in indepvar]
677 else:
678 return self.output(indepvar)
679 except (OverflowError, ValueError):
680 self.diagnostics.errors.append((indepvar, self.name + ": Overflow error in output"))
681 raise
682 except PyDSTool_BoundsError:
683 self.diagnostics.errors.append((indepvar, self.name + ": Bounds error in output"))
684 raise
685 elif checklevel in [1,2]:
686 if self.trajirange is None:
687 idep = self.indepdomain
688 else:
689
690 idep = self.trajirange
691 indepvar_ok = True
692
693
694 if isinstance(indepvar, _seq_types):
695 vectorizable = self._vectorizable
696 for d in indepvar:
697
698
699 try:
700 contresult = d in idep
701 except PyDSTool_UncertainValueError:
702 contresult = True
703
704
705 if d < idep[0]:
706 try:
707
708 dix = indepvar.index(d)
709 except AttributeError:
710
711 dix = indepvar.tolist().index(d)
712 indepvar[dix] = idep[0]
713 elif d > idep[1]:
714 try:
715
716 dix = indepvar.index(d)
717 except AttributeError:
718
719 dix = indepvar.tolist().index(d)
720 indepvar[dix] = idep[1]
721 if checklevel == 2:
722 self.diagnostics.warnings.append((d, None))
723 if not contresult:
724 indepvar_ok = False
725 break
726 else:
727 vectorizable = True
728 try:
729 indepvar_ok = indepvar in idep
730 except PyDSTool_UncertainValueError, errinfo:
731
732
733 if indepvar < idep[0]:
734 indepvar = idep[0]
735 elif indepvar > idep[1]:
736 indepvar = idep[1]
737 if checklevel == 2:
738 self.diagnostics.warnings.append((indepvar, None))
739
740
741 if not indepvar_ok:
742
743
744 if checklevel == 2:
745 self.diagnostics.errors.append((indepvar,
746 self.name + " : " + self.indepdomain._infostr(1)))
747 if vectorizable:
748 raise ValueError('Independent variable value(s) '
749 'out of range in Variable call')
750 else:
751 raise ValueError('Independent variable value '+\
752 str(indepvar) + ' out of '
753 'range in Variable call')
754 try:
755 if vectorizable:
756 depvar = self.output(indepvar)
757 else:
758 depvar = [self.output(ival) for ival in indepvar]
759 depvar_ok = True
760 except PyDSTool_BoundsError, errinfo:
761 depvar_ok = False
762
763 if depvar_ok:
764
765
766
767 if isinstance(depvar, (_seq_types, Pointset)):
768 if isinstance(depvar, Pointset):
769 dv = depvar.toarray()
770 else:
771 dv = depvar
772 for d in dv:
773
774
775 try:
776 contresult = d in self.depdomain
777 except PyDSTool_UncertainValueError, errinfo:
778 contresult = True
779 if checklevel == 2:
780
781
782 try:
783
784 depix = dv.index(d)
785 except AttributeError:
786
787 depix = dv.tolist().index(d)
788 self.diagnostics.warnings.append((indepvar[depix], errinfo.value))
789 if not isfinite(d):
790
791
792
793 raise PyDSTool_BoundsError("Return value was not finite/defined (%s)"%str(d))
794 if not contresult:
795 depvar_ok = False
796 break
797 elif depvar is None:
798
799
800
801
802 raise ValueError("Cannot compute a return value for "
803 "independent variable value "
804 + str(indepvar))
805 else:
806 if isinstance(depvar, Point):
807 dv = depvar[0]
808 else:
809 dv = depvar
810 try:
811 depvar_ok = dv in self.depdomain
812 except PyDSTool_UncertainValueError, errinfo:
813 if checklevel == 2:
814 self.diagnostics.warnings.append((indepvar, errinfo.varval))
815 if not isfinite(dv):
816
817
818
819 raise PyDSTool_BoundsError("Return value was not finite/defined (%s)"%str(dv))
820
821 if depvar_ok:
822 return dv
823 else:
824
825
826
827
828 if vectorizable:
829
830
831 raise PyDSTool_BoundsError('Computed value(s) %f outside'%dv + \
832 ' validity range in Variable call')
833 else:
834 raise PyDSTool_BoundsError('Computed value %f outside'%dv + \
835 ' validity range in Variable call')
836 else:
837
838 indepvar_ok = False
839 try:
840
841
842 if isinstance(indepvar, _seq_types):
843 vectorizable = self._vectorizable
844 indepvar_ok = all([i in self.indepdomain for i in \
845 indepvar])
846 else:
847 vectorizable = True
848 indepvar_ok = indepvar in self.indepdomain
849 except TypeError, e:
850 raise TypeError('Something messed up with the Variable '
851 'initialization: '+str(e))
852 else:
853 if not indepvar_ok:
854 raise ValueError('Independent variable '+str(indepvar)+\
855 ' out of range in Variable call')
856
857
858
859
860 try:
861 if vectorizable:
862 depvar = self.output(indepvar)
863 depvar_ok = depvar in self.depdomain
864 else:
865 depvar = [self.output(ival) for ival in indepvar]
866 depvar_ok = all([d in self.depdomain for d in \
867 depvar])
868 except PyDSTool_BoundsError, e:
869 raise ValueError("Cannot compute a return value for "
870 "this independent variable value: "
871 + str(e))
872 except PyDSTool_TypeError:
873 if not self.defined:
874 print "Variable '%s' not fully defined."%self.name
875 return None
876 else:
877 raise
878 else:
879 if depvar_ok:
880 return depvar
881 else:
882 if vectorizable:
883 raise PyDSTool_BoundsError('Computed value(s) '
884 'outside validity range in Variable call')
885 else:
886 raise PyDSTool_BoundsError('Computed value '+str(depvar)+\
887 'outside validity range in Variable call')
888
889
892
893
894 __str__ = __repr__
895
896
898 if verbose == 0:
899 return "Variable "+self.coordname+"("+self.indepvarname+")"
900 else:
901 try:
902 if isinputcts(self):
903 ipstr = "continuous"
904 else:
905 ipstr = "discrete"
906 except ValueError:
907 ipstr = "not defined"
908 outputStr = "Variable:\n Independent variable '" \
909 + self.indepvarname + "' [" + ipstr + "]\n"
910 try:
911 if isoutputcts(self):
912 opstr = "continuous"
913 else:
914 opstr = "discrete"
915 except ValueError:
916 opstr = "not defined"
917 outputStr += " defined in domain " + str(self.indepdomain)
918 if verbose == 2:
919 if self.trajirange is None:
920 outputStr += "\n ranges not known for this trajectory"
921 else:
922 outputStr += "\n trajectory ranges "+str(self.trajirange)
923 outputStr += "\nDependent variable '" + self.coordname + \
924 "' [" + opstr + "]\n defined in domain "
925 if not isinstance(self.depdomain, Interval):
926 outputStr += _num_type2name[self.coordtype]+": "
927 outputStr += str(self.depdomain)
928 if verbose == 2:
929 if self.trajdrange is None:
930 outputStr += "\n ranges not known for this trajectory"
931 else:
932 outputStr += "\n trajectory ranges "+str(self.trajdrange)
933 return outputStr
934
935
936 - def info(self, verboselevel=1):
938
939
941 pickledself = pickle.dumps(self)
942 return pickle.loads(pickledself)
943
944
946 pickledself = pickle.dumps(self)
947 return pickle.loads(pickledself)
948
949
951 d = copy.copy(self.__dict__)
952
953 d['indepvartype'] = _num_type2name[self.indepvartype]
954 d['coordtype'] = _num_type2name[self.coordtype]
955 if 'funcspec' in self._funcreg:
956
957
958
959
960 del d['output']
961 for fname, finfo in self._funcreg.iteritems():
962 if finfo[0] == 'self':
963 try:
964 del d[fname]
965 except KeyError:
966 pass
967
968
969
970 return d
971
972
974 self.__dict__.update(state)
975
976
977 self.indepvartype = _num_name2type[self.indepvartype]
978 self.coordtype = _num_name2type[self.coordtype]
979
980 for fname, finfo in self._funcreg.iteritems():
981 if finfo[0] == 'self' and not hasattr(eval(finfo[0]), fname):
982
983 setattr(eval(finfo[0]), fname, finfo[1])
984 if 'funcspec' in self._funcreg:
985
986 funcspec = self._funcreg['funcspec'][1]
987 outputdata = self._funcreg['outputdata'][1]
988 if hasattr(self, '_var_namemap'):
989 var_namemap = self._var_namemap
990 else:
991 var_namemap = None
992 if hasattr(self, 'initialconditions'):
993 ics = copy.copy(self.initialconditions)
994 else:
995 ics = None
996 if hasattr(self, '_refvars'):
997 if self._refvars is not None and self._refvars != []:
998 refvars = [copy.copy(v) for v in self._refvars]
999 else:
1000 refvars = None
1001 else:
1002 refvars = None
1003
1004 self.setOutput(outputdata, funcspec,
1005 self.globalt0, var_namemap, ics, refvars)
1006
1007
1009
1010
1011
1012
1013 for fname, finfo in self._funcreg.iteritems():
1014
1015 if finfo[0] is None:
1016
1017 continue
1018 elif fname == '_impfn':
1019 exec_str = 'del Variable.' + finfo[0]
1020 try:
1021 exec exec_str
1022 except AttributeError:
1023
1024
1025
1026 pass
1027 elif fname is 'funcspec':
1028
1029
1030 pass
1031 elif fname is 'outputdata':
1032
1033
1034 pass
1035 elif hasattr(eval(finfo[0]), fname):
1036 exec_str = 'del '+ finfo[0] + '.' + fname
1037 try:
1038 exec exec_str
1039 except RuntimeError:
1040
1041
1042 pass
1043 if hasattr(self, '_refvars'):
1044 if self._refvars is not None and self._refvars != []:
1045 for v in self._refvars:
1046 v.__del__()
1047
1048
1050 """Mimics part of the API of a non-hybrid variable.
1051
1052 This is a somewhat ugly hack as it's implemented by using a whole
1053 HybridTrajectory object to extract individual variable values,
1054 rather than having extracted a sequence of Variable objects from
1055 a HT and stitching them back together as a single entity."""
1056 - def __init__(self, hybridtraj, coordname, indepdomain, abseps=None):
1057
1058 self._ht = hybridtraj
1059 self.name = 'Hybrid variable '+coordname
1060 self.outputdata = None
1061 self.defined = True
1062 self.indepvarname = 't'
1063 self.indepdomain = indepdomain
1064 self.indepvartype = float
1065 self.coordname = coordname
1066 self.depdomain = Interval(self.coordname, float,
1067 [-Inf, Inf], abseps=abseps)
1068 self.coordtype = float
1069 self.trajirange = None
1070 self.trajdrange = None
1071 self.diagnostics = Diagnostics()
1072
1073
1074 self.output = None
1075
1076 - def __call__(self, indepvar, checklevel=0):
1077 return self._ht(indepvar, self.coordname, checklevel=checklevel)
1078
1080 """Returns a Pointset of independent and dependent variable values,
1081 provided variable is based on a mesh (otherwise None is returned).
1082 """
1083 return self._ht.sample([self.coordname])
1084
1086 """Reveal underlying mesh as arrays, rather than Pointset as returned
1087 by getDataPoints method."""
1088 vs = self._ht.sample([self.coordname])
1089 return array([vs.indepvararray, vs.coordarray[0]])
1090
1092 return "Hybrid variable "+self.coordname
1093
1094 __str__ = __repr__
1095
1096 - def info(self, verboselevel=1):
1097 return "Hybrid variable "+self.coordname
1098
1099
1100
1103
1105 self.__dict__.update(state)
1106
1110
1111
1112
1114 """One-dimensional function wrapper."""
1115
1116 - def __init__(self, fn, datapoints=None, numtypes=(float64,float64),
1117 abseps=None):
1118 assert isinstance(fn, types.FunctionType) or \
1119 isinstance(fn, types.BuiltinFunctionType), \
1120 ("fn argument must be a regular Python function")
1121 self.fn = fn
1122
1123
1124 if datapoints is None:
1125 datapoints = (Interval('indepvardom', numtypes[0], [-Inf, Inf],
1126 abseps=abseps),
1127 Interval('depvardom', numtypes[1], [-Inf, Inf],
1128 abseps=abseps))
1129 try:
1130 self.datapoints = (datapoints[0], datapoints[1])
1131 except TypeError:
1132 raise TypeError("datapoints argument must be a 2-tuple or list "
1133 "of 2-tuples or lists")
1134 try:
1135 self.types = (numtypes[0], numtypes[1])
1136 except TypeError:
1137 raise TypeError("numtypes argument must be a 2-tuple or list "
1138 "of 2-tuples or lists")
1139
1140
1142 if isinstance(arg, _seq_types):
1143 try:
1144 return self.fn(arg)
1145 except:
1146 return array([self.fn(v) for v in arg])
1147 else:
1148 return self.fn(arg)
1149
1150
1157
1158
1164
1165
1166
1167
1168
1187
1188
1191
1192
1193
1194
1195
1196
1197
1198
1199
1200
1201
1203 assert isinstance(var, Variable), "Argument must be a Variable"
1204 if var.defined:
1205 if compareNumTypes(var.coordtype, float64):
1206 return isinstance(var.depdomain, Interval) and not \
1207 isinstance(var.output, Pointset)
1208 elif compareNumTypes(var.coordtype, int32):
1209 return False
1210 else:
1211 raise TypeError("Unsupported dependent variable type for Variable")
1212 else:
1213 raise ValueError("Variable is not fully defined")
1214
1215
1218
1219
1221 """Determine if variable is continuously defined on its input and
1222 output domains."""
1223 assert isinstance(var, Variable), "Argument must be a Variable"
1224 return isinputcts(var) and isoutputcts(var)
1225
1226
1228 """Determine if variable is discretely defined on its input and
1229 output domains."""
1230 return not (isinputcts(var) and isoutputcts(var))
1231