o Bac
Some Dynamics Balance of Production via
Optimization and Simulation within System
Dynamics Method
Elzbieta Kasperska, Elwira Mateja-Losa, Damian Stota
Institute of Mathematics
Silesian University of Technology
Kaszubska 23, Gliwice 44-100, Poland
e-mail: clakaspe@polsl.gliwice.pl, elwimat@polsl.gliwice.pl, damslota@polsl-gliwice.pl
Abstract
The purpose of this paper is to present the dynamic optimal balance of production. Two
different formulations of this problem discussed here on the examples presented by authors.
The constrained and unconstrained optimization are presented on the basis of many simula-
tion experiments carried out by authors.
Keywords: System Dynamics Method, Optimization during Simulation, Simulation du-
ring Optimization, Balance of Production.
1 Introduction
The problem of optimal dynamics balance of flows is almost completly new on the area of System
Dynamics method. The first attempts in this subject were undertaken by Kasperska (Kasperska
1990), then by Kasperska and Slota (Kasperska and Slota 2000) and the newest publication
by Kasperska, Mateja-Losa and Slota (Kasperska et al. 2000a, Kasperska et al. 2000b). These
attempts were conected with unconstrained optimization during simulation and simulation du-
ring optimization in sense of Coyle (Coyle 1996, Coyle 1998, Coyle 1999). Now authors extended
these formulations to constrained versions of the problem. Many interesting experiments were
taken and some conclusions are presented.
2 The constrained and unconstrained optimizations
of the dynamics balance of production
To approach this problem we have to return to the example of dynamics balance of production,
discussed in paper (Kasperska et al. 2000a, Kasperska et al. 2000b). In Figure 1 an example of the
production system, in view of extended Lukaszewicz symbols (Lukaszewicz 1975, Lukaszewicz
1976), is presented.
The element in double surrounding” was called ” optimal balance of a, 3, 7”, and has more
important meaning here. Now we want to concentrate on these parts of the model which are
related to constrained and unconstrained optimization in sense of Coyle (Coyle 1996, Coyle
1999).
The parameters a, 3, 7 were optimized during simulation (and have lower and upper limits).
We assume that there are six cases of balances:
I) unconstrained balance of flows, a + 3 + 7 is optional (any), the formulation of the
objective function f,, is clear from the figure 1 (it consists of three elements with three
weighting factors);
Pj
f Y production
of item P1
production
of item P2
production
of item P3
source of
raw material
(input)
optimal
balance
of a, B,y
level of material
during transformation
ULPI,2,3,
—o ”
"FLB)
=} CPLA,
oo
st balance (SFFCB
the cos
Figure 1. Optimal dynamics balance of production.
II) constrained balance of flows, a + 3+ 7 = 1, this condition denotes full accordance
with the actual production of three items with optimal value of their production. Addition
of penalty function was required what resulted in the occurence of discrepancies from the
condition (a+3+7 = 1). When the condition is not fulfiled (a+8+y < lora+$+y7> 1)
the massive penalty factor named kara is added to the value of the base function fp.
Technically speaking it has a form of:
penalty = kara * max(0, abs(a+ 6+ —- 1)):
III) constrained balance of flows, a + 3+ 7 < 1. The condition of balance is now “not
sharp”. The difference between case II a this one is that now penalty function has the
form of:
penalty = kara * max (0, at+B+7- 1).
Only when a+ +7 > 1 the kara is added to the base function fo;
IV
constrained balance of flows, a+ 3+ 7 < 1. The condition of balance is now ”sharp”.
The logistic alternative function clip (see Coyle 1977, Coyle 1996, Forrester 1961, Forrester
1972) is now added:
penalty = kara « clip(1, Oa+B+7, 1).
The interesting point of view is: how to interpret such a condition a + 6+ 7 < 1. We
think that in such a case the actual values of possible production are bigger than those
required by optimal balance of production. In the affect it should cause the limitation of
actual production;
V) constrained balance of flows, a + 3+ 7 > 1. The condition of balance is now "not
sharp”. Similarly to IV, the logistic function clip is added:
penalty = kara clip(0,1,a+ 6+ 7,1).
We can interpret such a condition in the following way: the possibilities of production of
three items are too small when compared to those suggested by optimal balance (a+3+7 >
1) and this should cause the development of these possibilities;
VI) constrained balance of flows, a + 3+ > 1. This is a “sharp” version of type V.
Similarly to IV and V, we added the logic function clip:
penalty = kara * clip(1, O,la+ 6+ 7).
Different cases of balances (I-VI) require different base conditions of parameters a, 3, (see
tables 1-6 in the next section of this paper).
Many interesting experiments were undertaken by authors using COSMIC and COSMOS
(Cosmic 1994, Coyle 1996). We have to stress that we are on the beginning of the way to better
experimentation with COSMIC and COSMOS. Lots of attention should be paid to chosing the
factor kara and the base value of parameters.
Now we want to apply the constrained optimization to our idea presented in Lozanna (Ka-
sperska et al. 2000a, Kasperska et al. 2000b), named ” optimization during simulation”. In that
article we promised (in conclusions) that it would be interesting to investigate the so called
»pseudosolution” of differencies Mx —b (see Legras 1974) at the condition x; > 0, i 3. In
such a case of optimal balance of production we have to solve the system of equations, which is
created from the balance of the value of three properties of flows: mass balance ("rate of flow”
in Forrester sense), cost balance and personal balance. The idea of constrained optimization
(x; > 0) required the extension of matrix M to the form of:
a 1 1
ucp, ucpz UCp3
ulp, ulp2 ulps
M= 1 0 o |?
o 1 0
0 o 1
where: ucp; — unit cost of the production of product p;, i = 1,2,3, ulp; — unit labour of the
production of product p;, i = 1,2,3. The matrix b should be extended to the form:
b= (ro, tep, tle, ba, bs.bs)
where: rz — output rate of raw material (production), tcp — total cost expenditure of production,
tle — total labour expenditure of production, b;, i = 4,5,6 — the number of large value (chosen
experimentally). The solution of equation:
M-rx=b
has taken the form of (Legras 1974):
x= (M™.M)*-M?-b.
To computer programing of this we used DYNAMO for Windows (Dynamo 1994). Authors’
ideas were supported by their mathematical background and possibilities of language DYNAMO.
The readers can see this in Appendix A.
The comparison of the results of both authors’ investigations: optimization during simula-
tion (in sense of solving system of balance equation during simulation) and simulation during
optimization (in sense of Coyle), we can see discussing figures in the next section of the paper.
3 The results of simulation
First we want to present some of the results of the experiments that we achieved from simulation
during optimization (using COSMIC and COSMOS).
Tabel 1 contains the results of the first experiment type I of optimization (described in
section 2 of our paper). We can see that all of the three parameters: a, 3, y are assuming their
upper limits and find the value of function f,, which constitutes 34 % of its initial value.
Table 1.
Parameter | Final | Orginal | Lower Upper
value | value limit limit
a 1.0 0.50 0 ai
B 1.0 0.25 0 1
i 1.0 0.25 0 1
Initial value of objectiv function f,, | 0.124-10™
Final value of objectiv function f,, 0.420 - 10°
Final value of sf ftb 10.3 - 10‘
Final value of sf fcp 10.2 - 10
Final value of sf flp 80.5 - 10°
Table 2 contains the results of the experiment number 2 type I of optimization, which differ
from the assumption of the experiment number 1 about the values of parameters tcp and tle
(tep differ from 7000 to new value 3500 and tle from 1800 to 1000). Parameter a assumes its
lower imit, and final value of function f,, has 19 % of its initial value.
Table 2.
Parameter | Final | Orginal | Lower Upper
value | value limit limit
a 0.000 0.50 0 1
B 0.855 0.25 0 1
i 0.924 0.25 0 1
Initial value of objectiv function f,» | 0.203 - 10°
Final value of objectiv function f,4 0.383 - 10
Final value of sf ftb 15.7 - 10%
Final value of sf fcp 76.0 - 10°
Final value of sf flp 18.9 - 10°
Table 3 contains the results of the experiment number 3 of type IJ of optimization (a+8+ 7 =
1). Parameters a and ¥ assume their lower limits and parametr ( assumes its upper limit. The
value of objective function f,, is 61 % of its initial value. The characteristic result is that
sf ftb=0.
Table 4 contains the results of the experiment number 4 of type III of optimization (a+ 8+
and ¥ is considerably smaller than one, and the
7 <1). The sum of values of parameters a, {
value of objective function f,, is 84 % of its initial value.
Table 5 contains the results of the experiment number 5 of type IV of optimization (a +
+7 <1). The assumed values of a, 3 and ¥ are such that their sum is smaller than one. The
value of objective function f,, is 92 % of its initial value. Such results show that chosen initial
values were very close to the optimum values.
Table 3.
Parameter | Final | Orginal | Lower Upper
value | value limit limit
a 0.0 0.50 0 ai
B 1.0 0.25 0 1
7 0.0 0.25 0 1
Initial value of objectiv function f,, | 0.186 - 10%
Final value of objectiv function f,, | 0.113 - 10
Final value of sf ftb 0.0
Final value of sf fcp 20.89 - 10'
Final value of sf flp 71.74 - 10°
Table 4.
Parameter | Final | Orginal | Lower Upper
value | value limit limit
a 0.429 0.50 0 1
B 0.559 0.25 0 1
y 0.003 0.25 0 1
Initial value of objectiv function f,, | 0.124 -10™
Final value of objectiv function f,, 0.105 - 10
Final value of sf ftb 18.34 - 10
Final value of sf fcp 29.36 - 10
Final value of sf flp 73.23 - 10°
Table 5.
Parameter | Final | Orginal | Lower Upper
value | value limit limit
a 0.365 0.30 0 1
B 0.317 0.25 0 1
9 0.317 0.25 0 1
tiv function fy | 0.194 - 10
Final value of objectiv function f,, 0.176 - 10%
Initial value of obje
Final value of sf ftb 1.0314
Final value of sf fep 32.87 - 108
Final value of sf flip 97.95 - 10°
Table 6 contains the results of the experiment number 6 of type VI of optimization (a+ 6+
7 > 1). The parameters a, 3 and 7 assumed their upper limit values. The value of objective
function fo is ¢ % of its initial value.
Some of these results are presented in an illustrative form on figures 2-10.
Now we want to present some of the results of the experiments that we achieved from
optimization during simulation (using DYNAMO for Windows). Some values of parameters ucp1
(see Section 2 of this paper), tle, tep and values of bj, i = 4,5,6 vary in respective experiments
numbered 1-5 (see program in Appendix A). Th em appeared to be very sensitive to these
can see in figures 11-13.
values, what reader:
Table 6.
Parameter | Final | Orginal | Lower Upper
value | value limit limit
a 1.0 1.00 0 1
B 1.0 0.25 0 1
1.0 0.25, 0 1
Tnitial value of objectiv function f,, | 0.168 - 10™
Final value of objectiv function f,, 0.577 - 10°
Final value of sf ftb 10.32 - 10‘
Final value of sf fcp 10.22 - 10
Final value of sf flp 80.45 - 108
a 10.0 20.0 30.0 40.0 50.0, 60.0, 70.0 80,0 30.0 400.
TINE, <HEEK) TIME WITHIN SIMULATION FENINUH HAKTL
UNIT) LEVEL OF RAM MATERIAL DURING TRANSFORMATION 0.0
<UNIT/HEEK) OUTPUT RATE TO LMT 0
<UNIT/HEEK) INPUT RATE TO SFFTB 0.0
<$/MEEK) INPUT RATE TO SFFCB 5542.3
<NEN/MEEK) INPUT RATE TO SFFLB 666.3
Figure 2. The results of unconstrained optimization: a + 3+ = optional (the dynamics of the
characteristic of some variables of the production system).
0 10.0 20.0 30.0 40.0 50.0 60.0
TINE (EEK) TIME HITHIN SIMULATION
(1 OBJUCTION FUNCTION
CUNIT##2/HEEK> SUM FUNCTION OF FITTING TOTAL BALANCE
($9ee2/HEEK) SUM FUNCTION OF FITTING COST BALANCE
MEN®#2/HEEK) SUM FUNCTION OF FITTING PERSON BALANCE
Figure 3. The results of unconstrained optimization: a + 3+ = optional (the dynamics of the
characteristic of the objective function and its elements).
70.0
80.0 30.0 100.¢
PINInun MAXIMUM
O.00E+08 —1,19E+03
.00E+03 © 0,00E+09
O.00E+03 © 3,48E+09
0.00E+03 —0,09E+03
"ble 10.0 20.0 30.0 40.0 50.0 60.0
TIME WEEK) TIME WITHIN SIMULATION
— RPI UNIT/HEEK) RATE OF PRODUCTION OF FIRST ITEM (P1)
[ovr RP2 CUNIT/MEEK) RATE OF PRODUCTION OF SECOND ITEM <P2>
RP3 CUNIT/HEEK) RATE OF PRODUCTION OF THIRD ITEM <P3)
Figure 4. The results of unconstrained optimization: a + 3+ = optional (the dynamics of the
characteristic of the production rates for three item:
70.0
0.0
0.0
0.0
80.0 30.0
MINIMUM” MAXIMUM
100.0
323.9
162.0
162.0
TIME <HEEK) TIME HITHIN SIMULATION
UNIT) LEVEL OF RAH MATERIAL DURING TRANSFORMATION
CUNIT/HEEK) OUTPUT RATE TO LNT
CUNIT/HEEK) INPUT RATE TO SFFTB
<$/MEEK) INPUT RATE TO SFFCB
<MEN/MEEK) INPUT RATE TO SFFLE
Figure 5. The results of constrained optimization: a + 3 + = 1 (the dynamics of the charac-
teristic of some variables of the production system).
10.0 20.0 30.0 40.0 50.0 60.0 70.0 80.0 30.0 00.0
Fitncrun *° never
0.0
0.0
0.0
5542.3
666.3
1288.7
647.8
0.0
7000.0
1800.0
20.0 30.0 40.0 50.0 60.0 70.0 Er
TIME, WEEK) TINE WITHIN SIMULATION
(1) OBJUCTION FUNCTION
UNIT##2/HEEK) SUM FUNCTION OF FITTING TOTAL BALANCE
(ee2/HEEK) SUM FUNCTION OF FITTING COST BALANCE
<MENe#2/HEEK) SUM FUNCTION OF FITTING PERSON BALANCE
Figure 6. The results of constrained optimization: a + 3 + = 1 (the dynamics of the charac-
teristic of the objective function and its elements).
10.0,
MOnInunt
0.006403
0.00£+03
0.00£+03
0,00£+03
© ante
1.786+09
0.006403
3.48E+03
0.096403
100 200 000 40.0 600 «20 700 00.0 90.0, 1900
TIME (HEEK) TIME MITHIN SIMULATION Te eH
CUNIT/HEEK) RATE OF PRODUCTION OF FIRST ITEM <P1) o.0 323.8
CUNIT/HEEK) RATE OF PRODUCTION OF SECOND ITEM (P2) a0 162.0
CUNIT/HEEK) RATE OF PRODUCTION OF THIRO ITEM <P3> a.0 162.0
Figure 7. The results of constrained optimization: a + 3 + = 1 (the dynamics of the charac-
teristic of the production rates for three items).
4200+
20.0 30.0 40.0 50.0 60.0 70.0
80.0 30.0 400.C
TIME, MEEK) TIME ITHIN SIMULATION TaN hurr
UNIT) LEVEL OF RAH MATERIAL DURING TRANSFORMATION 0.0 1295.7
CUNIT/HEEK) QUTPUT RATE TO LMT 0.0 647.8
CUNIT/MEEK) INPUT RATE TO SFFTB 0.0 123.6
($AMEEK) INPUT RATE TO SFFCB 5671.9 7000.0
(NEN/HEEK) INPUT RATE TO SFFLE 925.4 1800.0
Figure 8. The results of constrained optimization: a + 3 + < 1 (the dynamics of the charac-
teristic of some variables of the production system).
10.0 20.0 30.0 40.0 50.0 60.0 70.0
20.0 30.0 400.0
TINE CHEEK) TINE MITHIN: SIMULATION PERERA
(1) OBJUCTION FUNCTION 0.006+03 —1,.86E+08
CUNIT#82/AIEEIO SUM FUNCTION OF FITTING TOTAL BALANCE 0.006+03 — 0,.006+08
SFFCP (Bex2/EEK) SUN FUNCTION OF FITTINS COST BALANCE O.00E+03 3. 58E+08
“5 SFFLP<CMENss2/HEEK) SUM FUNCTION OF FITTING PERSON BALANCE O.00E+03 — 0.13E+08
Figure 9. The results of constrained optimization: a + 3+ < 1 (the dynamics of the charac-
teristic of the objective function and its elements).
20
8 10.0 20.0 30.0 40.0 50.0 60.0 70.0 80.0 30.0 100.0
TIME CHEEK) TIME HITHIN SIMULATION POCA HERE
“ppt CUNIT/EEK> RATE OF PRODUCTION OF FIRST ITEM <P1) 0.0 194.4
[ie RP2 — CUNIT/HEEK RATE OF PRODUCTION OF SECOND ITEM <P2> a0 162.0
t==- RP3 CUNIT/HEEK) RATE OF PRODUCTION OF THIRD ITEM (P3) 0.0 162.0
Figure 10. The results of constrained optimization: a + 3+ <1 (the dynamics of the charac-
teristic of the production rates for three items).
menRate 11 men ate 1p em ate 1p3
Song eae 2 es Rate rp2 wee Le vel Imt
1500.
1000.
0 10 20 30 40 50 60 70 80 30 104
TIME
Figure 11. The dynamics of the characteristic of some variables of the model by optimization
during simulation (ucp1 = 10, tep = 4000, tle = 990, ba = bs = bg = 300).
oaeRatorl Rate 1p! men Rate 1p3
2000 mate 12 Rate 1p2 Le vel Imt
1500
1000
500,
a
0 10 20 30 40 50 60 70 80 30 104
TIME
Figure 12. The dynamics of the characteristic of some variables of the model by optimization
during simulation (ucp; = 10, tcp = 3500, tle = 1000, b4 = bs = bg = 50).
meeeRate 11 meen ate 1p mms ate 1p3
2009 Sate 2 es Rate rp2 wwe Le vel Imt
1500.
1000.
500.
0 10 20 30 40 50 60 70 80 30 104
TIME
Figure 13. The dynamics of the characteristic of some variables of the model by optimization
during simulation (ucp; = 10, tep = 4000, tle = 990, b4 = 100, bs = bg = 50).
4 Conclusion
After modelling and simulating some optimal balance of production we have come up to the
following conclusions:
a) The optimization during simulation allows to achieve the optimal pseudosolutions of the
underdetermined system of equations (in sense of Legras (Legras 1974)). The solutions
keep limitation x; > 0, i = 1,2,3, which was required by the physical sense of flows. The
future experiments will go towards more precise selection of the value of ” large” variables
(see section 2) in matrix b, to minimize the norm of the discrepancies (Ax; —b;), i = 1, 2,3.
b
The simulation during optimization (in sense of Coyle) allows to achieve the optimal
solution (the value of parameters a, 3, y) and the objective function f,». The future
experiments should extend the set of optimizing parameters, for example the parameters:
tle, tep, ucp;, i = 1,2,3. Much more attention should be paid to selecting the base value
of parameters a, 3, 7 in different kinds of constrained optimization.
The authors are on the beginning of the way to better experimentation with COSMIC and
COSMOS. We are planning to extend the optimal balancing of flows to more value of dimension
and more sophisticated form of the objective function.
We want to thank to Prof. R. G. Coyle for the inspiration from his books and articles and
for his COSMOS and COSMIC which are the fruitful tools for working with System Dynamics
models.
References
COSMIC and COSMOS user manuals. 1994. The COSMIC Holding Co.: Shrivenham.
Coyle, R. G. 1977. Management System Dynamics. John Wiley & Sons: New York.
Coyle, R. G. 1978. System Dynamics — The state of the art. Dynamica 5: 3-23.
Coyle, R. G., Wolsterholm, E. P. 1980. Modelling discrete events in System Dynamics model.
A case study. Dynamica 6: 21-27.
Coyle, R. G. 1996. System Dynamics modelling. A practical approach. Chapman & Hall: Lon-
don.
Coyle, R. G. 1998. The practice of System Dynamics: milestones, lessons and ideas from 30 years
experience. System Dynamics Review 14: 343-365.
Coyle, R. G. 1999. Simulation by repeated optimization. Journal of the Operational Research
Society 50: 429-438.
Forrester, J. W. 1961. Industrial Dynamics, MIT Press: Massachusetts.
Forrester, J. W. 1972. Principles of Systems, Cambridge Press: Massachusetts.
Forrester, J. W. 1975. Collected papers of Jay W. Forrester. Cambridge Press: Massachusetts.
Kasperska, E. 1990. Methods of simulation of the investigation into supporting planning and
organization in industry with continuous processes. Ph.D. thesis. Polish Academy of Science:
Warsaw (in Polish).
Kasperska, E., Slota, D. 2000. Mathematical Method of Management in the Concept of System
Dynamics. Silesian Technical University: Gliwice (in Polish).
Kasperska, E., Mateja-Losa, E., Slota, D. 2000. Some extension of System Dynamics method
theoretical aspects. In Proc. 16th IMACS World Congress. Deville M. and Owens R. (eds).
IMACS: Lausanne; 718-10, 1-6.
Kasperska, E., Mateja-Losa, E., Slota, D. 2000. Some extension of System Dynamics method
practical aspects. In Proc. 16th IMACS World Congress. Deville M. and Owens R. (eds). IMACS:
Lausanne; 718-11, 1-6.
Legras, J. 1974. Praktyczne metody analizy matematycznej. WNT: Warszawa (translation from
French: Methodes et Technique De l’Analyse Numerique. Dunod: Paris 1971).
Lukaszewicz, R. 1975. Management System Dynamics. PWN: Warsaw (in Polish).
Lukaszewicz, R. 1976. The direct form of structure models within System Dynamics. Dynami-
ca 2: 36-43.
Professional DYNAMO 4.0 for WINDOWS. Reference manual. 1994. Pugh-Roberts Associates:
Cambridge.
Appendix A. Program in DYNAMO
* Balance of production of 3 items
note
note level of raw material during transformation
k=1mt .jtdt*(r1.jk-r2. jk)
note
note input rate to lmt (r1)
note
vr ri.kl=g1l*input.k
c gi=5
note
note input - source of raw material
note
a input .k=pot+p1*sin((6.28*time.k)/perd)
c po=100
c pl=30
c perd=52
note
note output rate from lmt (r2)
note
xv r2.kl=lmt.k/t1
ec ti=2
note
note unit cost of production of 1-st item (ucp1)
note unit cost of production of 2-st item (ucp2)
note unit cost of production of 3-st item (ucp3)
note
c¢ ucpi=10
c ucp2=5
© ucp3=2
note
note unit labour of production of 1-st item (ulp1)
note unit labour of production of 2-st item (ulp2)
note unit labour of production of 3-st item (ulp3)
note
ce ulpi=2
c ulp2=2
c ulp3=1
note
note total cost of expediture of production (tcp)
note total labour of expediture of production (tlc)
note
c tcp=4000
c tle=990
note
note vector b
note
a bi.k=r2.k1
a b2.k=tcp
a b3.k=tle
a b4.k=100
a b5.k=50
a b6.k=50
note
note vector at.b
note
a bb1.k=b1.k+b4.k+b2.k*ucp1+b3.k*ulp1
a bb2.k=b1.k+b5.k+b2.k*ucp2+b3.k*ulp2
a bb3.k=b1.k+b6.k+b2.k*ucp3+b3.k*ulp3
note
note matrix c=at.a
note
c11.k=2+ucp1*ucpi+ulpi*ulpt
c12.k=1+ucp1*ucp2+ulp1*ulp2
c13.k=1+ucp1*ucp3+ulp1*ulp3
c21.k=1+ucp1*ucp2+ulp1*ulp2
c22.k=2+ucp2*ucp2+ulp2*ulp2
c23.k=1+ucp2*ucp3+ulp2*ulp3
c31.k=1+ucp1*ucp3+ulp1*ulp3
c32.k=1+ucp2*ucp3+ulp2*ulp3
c33.k=2+ucp3*ucp3+ulp3*ulp3
note
sep ep pp ppp pw
note determinant of matrix c=at.a
note
a detc.k=-c13.k*c22.k*c31.k+c12.k*c23.k*c31.k+c13.k*c21.k*c32.k 7
~ci1.k*c23.k*c32.k-c12.k*c21.k*c33.k+c11.k*c22.k*c33.k
note
note matrix d=Det[c]*Inverse[c]
note
dii1.k=c22.k*c33.k-c23.k*c32.k
d12.k=c13.k*c32.k-c12.k*c33.k
d13.k=-c13.k*c22.k+c12.k*c23.k
d21.k=-c21.k*c33.k+c31.k*c23.k
d22.k=-c13.k*c31.k+c11.k*c33.k
d23.k=c13.k*c21.k-c11.k*c23.k
d31.k=-c22.k*c31.k+c21.k*c32.k
d32.k=c12.k*c31.k-c11.k*c32.k
d33.k=-c12.k*c21.k+c11.k*c22.k
note
note rate of production of 1-st item (rp1)
note rate of production of 2-st item (rp2)
note rate of production of 3-st item (rp3)
note
a rpi.k=(bb1.k*d11.k+bb2.k*d12.k+bb3.k*d13.k)/detc.k
a rp2.k=(bb1.k*d21.k+bb2.k*d22.k+bb3.k*d23.k)/detc.k
a rp3.k=(bb1.k*d31.k+bb2.k*d32.k+bb3.k*d33.k)/detc.k
note
sep ep ppp pp
note
note
bli.k=(b1.k-rp1.k-rp2.k-rp3.k)*(b1.k-rp1.k-rp2.k-rp3.k)
b12.k=(ucp1*rp1.k+ucp2*rp2.k+ucp3*rp3.k-b2.k)
*(ucpi*rp1.k+ucp2*rp2.k+ucp3*rp3.k-b2.k)
b13.k=(ulp1*rp1.k+ulp2*rp2.k+ulp3*rp3.k-b3.k)
*(ulpi*rp1.k+ulp2*rp2.k+ulp3*rp3.k-b3.k)
b1l4.k=(rp1.k-b4.k) *(rp1.k-b4.k)
b15.k=(rp2.k-b5.k) *(rp2.k-b5.k)
b16.k=(rp3.k-b6.k) * (rp3.k-b6.k)
functio.k=bl1.k+b12.k+b13.k+b14.k+b15.k+b16.k
note
note parameters of simulation
note
spec length=104/dt=1/savper=1
save rp1,rp2,rp3,r1,r2,lmt,detc,functio
pppek px p p
Appendix B. Program in COSMIC
note Balance of production of 3 items
1mt.k=Imt. j+dt*(r1.jk-r2.jk)
1mt=0
gi=5
po=100
p1=30
perd=52
t1=2
ri.kl=input.k*g1
tcp=7000
ucp1=10
ucp2=5
ucp3=2
alfa=0.5
beta=0.25
gamma=0.25
tle=1800
ulp1=2
ulp2=2
ulp3=1
wi=10
w2=2
w3=2
r2.kl=lmt.k/t1
input .k=potp1*SIN(6.283*time.k/perd)
fob. k=wi*sfftb.k+w2*sffcp.k+w3*sfflp.k
rsi.kl=(1-alfa-beta-gamma) *ar2.k
ar2.k=Imt.k/t1
rpi.kl=alfa*ar2.k
rp2.kl=beta*ar2.k
rp3.kl=gamma*ar2.k
rs2.kl=tcp-rp1.kl*ucp1-rp2.kl*ucp2-rp3.kl*ucp3
rs3.kl=tle-rp1.kl*ulp1-rp2.k1*ulp2-rp3.k1*ulp3
sfftb.k=sfftb.j+dt*((rs1.jk)**INT(2))
PRHHHHPHPHMHAAAAAANRAANXNAAAAAHAAAA ABB
n sfftb=wp1
1 sffcp.k=sffcp. j+dt*((rs2. jk) **INT(2))
n sffcp=wp2
1 sfflp.k=sfflp.j+dt*((rs3. jk) **INT(2))
n sfflp=wp3
c wp1=0
c wp2=0
c wp3=0
note
note output and control sector
note
ce dt=1
c length=104
c prtper=5
c pltper=5
print 1)rpi,rp2,rp3
plot rpi=1,rp2=2,rp3=3
note
note simulation experiments
note
run basic model of balance
note
note definitions of variables
note
alfa=(1) fraction of production ar2 (of item P1)
ar2=(unit/week) transformation of raw material (production)
beta=(1) fraction of production ar2 (of item P2)
dt=(week) solution interval
fob=(1) objuction function
gi=(1) gain between input and output rate r2
gamma=(1) fraction of production ar2 (of item P3)
input=(unit/week) source of raw material
length=(week) simulated period
lmt=(unit) level of raw material during transformation
pi=(1) parameter of amplitude of sinusoidal input
perd=(week) parameter, period of input
pltper=(week) graph plotting interval
po=(1) parameter step of input
prtper=(week) table printing interval
ri=(unit/week) input rate to Imt
r2=(unit/week) output rate to 1lmt
rpi=(unit/week) rate of production of first item (p1)
rp2=(unit/week) rate of production of second item (P2)
rp3=(unit/week) rate of production of third item (P3)
rsi=(unit/week) input rate to sfftb
rs2=(dolar/week) input rate to sffcb
rs3=(men/week) input rate to sfflb
sfficp=($**2/week) sum function of fitting cost balance
sfflp=(men**2/week) sum function of fitting person balance
sfftb=(unit**2/week) sum function of fitting total balance
ti=(week) average time of production
anaanaaAnAAAAAAAAAAAAAAAAARAARAAA
aeaanaAaAanAAAAAAAA
tcp=($/week) total cost expenditure of production (rate)
time=(week) time within simulation
tle=(men/week) total labour expenditure of production (rate)
ucp1=($/unit) unit cost of production of product P1
ucp2=($/unit) unit cost of production of product P2
ucp3=($/unit) unit cost of production of product P3
ulpi=(men/unit) unit labour of production of product P1
ulp2=(men/unit) unit labour of production of product P2
ulp3=(men/unit) unit labour of production of product P3
wi=(1) weighting factor for sfftb
w2=(1) weigting factor for sffcb
w3=(1) weighting factor for sfflb
wp1=(unit**2/week) initial condition for sfftb
wp2=($**2/week) initial condition for sffcb
wp3=(men**2/week) initial condition for sfflb
