0: #!/usr/bin/perl
1: use strict;
2: use warnings;
3:
4: # Perl solving a physics / electrodynamics problem involving
5: # symbolic mathematics, derivatives and complex numbers:
6:
7: use Math::Symbolic qw/:all/;
8: use Math::Complex;
9:
10: # Given the following simple circuit:
11: #
12: # /\/\/\ (R = resistor,
13: #  R L  L = solenoid,
14: #   U = alternating voltage)
15: # O ~ O
16: # U(t)
17: #
18: # Question: What's the current in this circuit?
19: #
20: # We'll need some physics before letting the computer do the
21: # math:
22: # Applying Kirchhoff's rules, one quickly ends up with the
23: # following differential equation for the current:
24: # (L * dI/dt) + (R * I) = U
25:
26: my $left = parse_from_string('L * total_derivative(I(t), t) + R *
+ I(t)');
27: my $right = parse_from_string('U(t)');
28:
29: <readmore>
30: # If we understand current and voltage to be complex functions,
31: # we'll be able to derive. ("'" denoting complex here)
32: # I'(t) = I'_max * e^(i*omega*t)
33: # U'(t) = U_max * e^(i*omega*t)
34: # (Please note that omega is the frequency of the alternating volt
+age.
35: # For example, the voltage from German outlets has a frequency of
+50Hz.)
36:
37: my $argument = parse_from_string('e^(i*omega*t)');
38: my $current = parse_from_string('I_max') * $argument;
39: my $voltage = parse_from_string('U_max') * $argument;
40:
41: # Putting it into the equation:
42: $left>implement( I => $current );
43: $right>implement( U => $voltage );
44:
45: $left = $left>apply_derivatives()>simplify();
46:
47: # Now, we can solve the equation to get a complex function for
48: # the current:
49:
50: $left /= $argument;
51: $right /= $argument;
52: my $quotient = parse_from_string('R + i*omega*L');
53: $left /= $quotient;
54: $right /= $quotient;
55:
56: # Now we have:
57: # $left = $right
58: # I_max(t) = U_max / (R + i*omega*L)
59: # But I_max(t) is still complex and so is the righthandside of t
+he
60: # equation!
61:
62: # Making the symbolic i a "literal" Math::Complex i
63: $right>implement(
64: e => Math::Symbolic::Constant>euler(),
65: i => Math::Symbolic::Constant>new(i), # Math::Complex magi
+c
66: );
67:
68: print <<'HERE';
69: Sample of complex maximum current with the following values:
70: U_max => 100
71: R => 10
72: L => 10
73: omega => 1
74: HERE
75:
76: print "Computed to: "
77: . $right>value(
78: U_max => 100,
79: R => 10,
80: L => 10,
81: omega => 1,
82: ),
83: "\n\n";
84:
85: # Now, we're dealing with alternating current and voltage.
86: # So let's make a generator that generates nice current
87: # functions of time!
88: # I(t) = Re(I_max(t)) * cos(omega*t  phase);
89:
90: # Usage: generate_current(U_Max, R, L, omega, phase)
91: sub generate_current {
92: my $current = $right>new(); # cloning
93:
94: $current *= parse_from_string('cos(omega*t  phase)');
95:
96: $current>implement(
97: U_max => $_[0],
98: R => $_[1],
99: L => $_[2],
100: omega => $_[3],
101: phase => $_[4],
102: );
103: $current = $current>simplify();
104: return sub { Re( $current>value( t => $_[0] ) ) };
105: }
106:
107: print "Sample current function with: 230V, 2Ohms, 0.1H, 50Hz, PI/
+4\n";
108: my $current_of_time = generate_current( 230, 2, 0.1, 50, PI / 4 )
+;
109:
110: print "The current at 0 seconds: " . $current_of_time>(0) . "\
+n";
111: print "The current at 0.1 seconds: " . $current_of_time>(0.1) .
+"\n";
112: print "The current at 1 second: " . $current_of_time>(1) . "\
+n";
113:
114: </readmore>
