Hindawi Publishing Corporation

Mathematical Problems in Engineering

Volume 2013, Article ID 781043, 7 pages
http://dx.doi.org/10.1155/2013/781043

Research Article
Long-Term Load Forecasting Based on a Time-Variant Ratio
Multiobjective Optimization Fuzzy Time Series Model
Xiaojuan Liu1,2 and Jianan Fang1
1
2

College of Information Science & Technology, Donghua University, Shanghai 201620, China
School of Mathematics and Physics, Shanghai University of Electric Power, Shanghai 201300, China

Correspondence should be addressed to Jianan Fang; jafang@dhu.edu.cn

Received 9 January 2013; Accepted 4 March 2013
Academic Editor: Engang Tian
Copyright 2013 X. Liu and J. Fang. This is an open access article distributed under the Creative Commons Attribution License,
which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Load forecasting problem is a complex nonlinear problem linked with economic and weather factors. Long-term load forecasting
provides useful information for maintenance scheduling, adequacy assessment, and limited energy resources for electrical power
systems. Fuzzy time series forecasting models can be used for long-term load forecasting. However, the interval length has been
chosen arbitrarily in the implementations of known fuzzy time series forecasting models, which has an important impact on the
performance of these models. In this paper, a time-variant ratio multiobjective optimization fuzzy time series model (TV-RMOP)
is proposed, and its performance is tested on the prediction of enrollment at the University of Alabama. Results clearly promote the
forecasting accuracy as compared to the conventional models. A genetic algorithm is used to search for the length of intervals based
on the training data while Pareto optimality theory provides the necessary conditions to identify an optimal one. The TV-RMOP
model is applied for the long-term load forecasting in Shanghai of China.

1. Introduction
Electric power load forecast has been a research topic for
many decades and the accuracy of load forecast is crucial
to electricity power industry due to its direct influence on
generating planning. During the past years, there have been
numerous attempts to improve the accuracy of load forecasting methods [1, 2]. The forecasting model and methods can be
divided into two categories: statistical methods and computational intelligence algorithms. The former method include
regression analysis and time series [3, 4]. The latter include
neural network, expert system, fuzzy logical, and so forth [5
7]. Now more and more scholar and researchers pay much
attention to uncertainty hybrid method [810].
There are three different types of electric load forecasting
depending on the time horizon and the operating decision
that needs to be made: short-, medium-, and long-term forecasting. In general, long term forecasting includes prediction
making for a range more than a year and is needed for power
supply and delivery system. Therefore, the weather factor is
not the main factor to affect the long-term load forecast, only
using historical load data can predict the 1 year ahead load.

Fuzzy time series [11, 12] prediction method is likely to get a

satisfactory result.
Fuzzy time series has been widely studied for recent years
for the aim of forecasting. The basic process of conventional
model is shown in Figure 1(a). Related studies mainly focus
on three factors: the partition of discourse, the content
of forecasting rules, and the methods of defuzzification,
all of which influence the prediction accuracy of forecasting
models. A number of fuzzy time series forecasting models
have been presented since Song and Chissom [13, 14] first
presented the concepts of fuzzy time series based on the
fuzzy set theory. However, the length of intervals has been
chosen arbitrarily in known models. This paper focuses
on the partition of discourse in TV-RMOP. We propose a
novel approach to determine the length of intervals based
on multiobjective optimization method in order to obtain
more accurate forecasts in fuzzy time series. The proposed
approach is based on two-parameter constrained multiobjective optimization. A genetic algorithm is used to search
for the two parameters based on the training data. The TVRMOP model is shown in Figure 1(b). Its performance is
tested on the prediction of enrollment at the University

Mathematical Problems in Engineering

Training data

Defining and partitioning the universe of

discourse arbitrarily
Defining fuzzy sets and fuzzifying time
series

Training data
Defining and partitioning the universe of
discourse with parameters
Defining fuzzy sets and fuzzifying time
series
Establishing fuzzy relationships

Establishing fuzzy relationships

Forecasting and defuzzifying

forecasting results

Forecasting and defuzzifying

forecasting results
Multiobjective optimization

Forecasting results

Testing data
(a) Conventional process

Forecasting results
Testing data
(b) Our process

Figure 1: Two distinct prediction processes of fuzzy time series forecasting.

of Alabama and applied for the long-term load forecasting

in Shanghai of China. Results investigated that TV-RMOP
model considerably improves the forecasting performance.
The rest of this paper is organized as follows. Section 2
introduces the methodology background briefly. The presented model is described, and the forecasting results on
the prediction of enrollment at the University of Alabama
are presented and discussed in Section 3. In Section 4, TVRMOP model is applied for the long-term load forecasting
in Shanghai of china. Finally, conclusions are presented in
Section 5.

2. Background
2.1. Fussy Time Series. The concepts of fuzzy time series are
described as follows [13, 14].
Let be the universe of discourse, where = {1 ,
2 , . . . , }. A fuzzy set of is defined as = (1 )/1 +
(2 )/2 + + ( )/ , where is the membership
function of the fuzzy set , : [0, 1]. is a generic
elements of fuzzy set . ( ) is the degree of belongingness
of to , ( ) [0, 1], and 1 .
Definition 1. Let () ( = . . . , 0, 1, 2, . . .) be the universe of
discourse and also a subset of real numbers. It is assumed that
() is defined on () and () is the collection of (), then
() is called a fuzzy time series on ().
Definition 2. It is assumed that () is a fuzzy time series and
() = ( 1) (, 1), where (, 1) is a fuzzy
relation and is an operator which is caused by ( 1).
The relationship between () and ( 1) can be denoted
by ( 1) () when () = ( 1) (, 1) is the
first-order fuzzy time series model of ().

Definition 3. Let () be a fuzzy time series. For any , () =

( 1) and () have only finite elements and therefore ()
is a time-invariant fuzzy time series; otherwise, it is a timevariant fuzzy time series.
2.2. Multiobjective Optimization and Pareto Points. It is common in decision-making problems to have a set of parameters
that play a dominant role in the process. It is also desirable
for the decision made to be the best possible with respect
to the available data. In the context of mathematics, the best
decision coincides with the optimization of a function of the
decision parameters, known as cost function or objective.
On the one hand, if the goal is to optimize a single cost
function, then the problem is classified as single objective
optimization. On the other hand, the process of simultaneously optimizing two or more conflicting objectives, subject
to certain constraints, is identified as multiobjective optimization [1719].
Such problems are called multiobjective optimization
problems and formulated as follows:
F (x) = [1 (x) , 2 (x) , (x)]
min
x
subject to

(x) 0,

= 0, 1, . . . , 1

(x) = 0,

= 0, 1, . . . , 2 ,

(1)

where denotes the number of objective functions; 1 , 2

identify the number of inequality and equality constraints,
respectively. x is an -dimensional vector in the solution
space X , and x = (1 , 2 , . . . , ) X .
Multiobjective optimization problems do not have a single global solution in contrast to the single objective problems. As a result, it is required to specify a number of criteriaconditions for which a set of solutions can be identified as

Mathematical Problems in Engineering

3
Start

Initialization of genetic
algorithm parameters
gaoptimset()
Assign new parameters
to TV-RMOP model

Calculate the fitness

function for the new
population

Generate a
new
population
using
reproduction,
mutation, and
Record the best results
crossover
mechanisms
No

Termination criterion?
Yes

Split samples into training

(Dtrain ) and test (Dtest ) sets

Defining and partitioning the

universe of discourse with two
parameters initial and ratio
Defining fuzzy sets and
fuzzifying time series

The intervals are calculated

by the following formulas:
upper0 = initial
for 1
lower = upper1
upper = (1 + ratio) upper0
interval = [lower , upper ]

Establishing fuzzy relationships

Forecasting and defuzzifying
forecasting results

Let ( 1) = , the defuzzified

forecast is equal to the midpoint
of the interval which corresponds
to

Multiobjective optimization
gamultiobj()

Perform forecasting for Dtest

End

Figure 2: Proposed time-variant ratio multiobjective optimization fuzzy time series model.

optimal. In general, the optimal solutions for F(x) can be

obtained based on the Pareto optimality theory. The set of
optimal solutions, which is known as the Pareto optimal set
X , are those solutions such that each of these cannot improve
any objective function in F(x) without the worsening of
some other objective functions at the same time. The optimal
solution set is defined as
X = {x X | does not exist x X, F (X) F (X )} . (2)
The corresponding points, which are known as the nondominated points, form the Pareto front in the objective space.
Since there exists different tradeoff solutions to minimize
F(x), one possible way to solve the multiobjective problem is
to select one solution manually from the Pareto optimal set.

Mean Absolute Error (MAE):

MAE () =

=1

Mean Absolute Percentage Error (MAPE):

MAPE () =

100

=1

Maximum Absolute Percentage Error (MAP):

MAP () = max (100
)

=1,...,
Theil Inequality Coefficient (THEIL):

3. Modeling and Verification

3.1. Proposed TV-RMOP Model. Time series is divided into
two parts. One part is used for training purpose while the
other part is used for testing. The proposed model was applied
to the training data to make estimation and then the predicted
values were found for the test data. In the training phase,
five measures are adopted to evaluate the model prediction
accuracy with the intention to formulate the multiobjective
problem. These measures are objective functions with parameters and defined as follows:

THEIL () =

RMSE
(
=1

( ) ) /

,
(3)

Root Mean Square Error (RMSE):

where is the actual value and is the predicted value at

time instant . is the number of predictions made for a
prespecified time interval. All measures should be minimized
in the optimization procedure.
As shown in Figure 2, the TV-RMOP model can be summarized as follows.

1
2
( )
=1

Step 1. Define the universe of discourse and intervals using

equations defined in the upper-right dashed box. Where
initial and ratio are two unknown parameters.

RMSE () =

Mathematical Problems in Engineering

Table 1: The enrollment data.

Step 2. Define fuzzy sets based on the universe of discourse

and fuzzify the historical data.

3.15
3.1
Objective 2

Enrollments
13055
13563
13867
14696
15460
15311
15603
15861
16807
16919
16388
15433
15497
15145
15163
15984
16859
18150
18970
19328
19337
18876

3.05
3
2.95
2.9
2.85
2.8
530

535

540

545

550 555
Objective 1

560

565

570

575

Figure 3: Pareto front of TV-RMOP model for enrollment.

Pareto front

5.96
5.94
5.92
Objective 2

Year
1971
1972
1973
1974
1975
1976
1977
1978
1979
1980
1981
1982
1983
1984
1985
1986
1987
1988
1989
1990
1991
1992

Pareto front

3.2

5.9
5.88
5.86
5.84
5.82

Step 3. Establish fuzzy relationships: .

Step 4. Forecast. Let (1) = ; the forecast is the midpoint
of the interval which corresponds to .

5.8
44

44.2

44.4

44.6
44.8
Objective 1

45

45.2

45.4

Step 5. Calculate RMSE, MAE, MAPE, MAP, and THEIL,

respectively. Each of them is a function with variable initial
and ratio.

Figure 4: Pareto front of TV-RMOP model for load forecasting.

In the training phase, the training samples (train ) are

used to calculate the above five functions. The determination
of variable initial and ratio is formulated as a multiobjective
optimization problem:

in the dashed box on the left in Figure 2. In this paper we use

a MATLAB function gamultiobj to optimize the initial and
ratio parameters.

min

(initial, ratio)
= [RMSE, MAE, MAPE, MAP, THEIL]

(4)

subject to initial (, ) , ratio (, ) .

In the test phase, an optimal solution is brought into the

model, and the forecast value is determined in the same way
as the training phase.
Genetic algorithm can be used to search for a solution of
multiobjective optimization problems. The process is shown

3.2. Verification of Proposed Model. In order to show what

we have achieved from the TV-RMOP model, we apply it
to the enrollment of University of Alabama and then have
compared the obtained results with those obtained from
Chens method [20]. Table 1 lists the enrollment from 1971
to 1992. The observations from 1971 to 1989 are used for the
training (train ), while the observations from 1990 to 1992 are
used for test.
The minimal value of train is 13055, then [12000, 13000]
is taken as constrains of variable initial. (0, 0.2] is taken as
constrains of the other variable ratio, because when the ratio
is greater than 0.2, the length of interval will be too large, and
the number of intervals will decrease. It results in diminishing

Mathematical Problems in Engineering

5
Table 2: Pareto solutions of enrollment.

Initial
12266.7
12230.4
12135.5
12846.0
12266.7

Ratio
0.0621
0.0530
0.0455
0.0676
0.0620

RMSE
571.9
573.9
532.7
568.0
551.5

MAE
484.4
486.0
451.8
501.2
459.0

MAPE
2.9745
2.9983
2.8276
3.1217
2.8142

MAP
5.8418
5.6954
6.6546
5.2926
5.8748

THEIL
0.0359
0.0360
0.0334
0.0357
0.0346

Table 3: The forecasting results for test .

Year
1990
1991
1992

Actual enrollment
19328
19337
18876

[15]
18685
19138
19176

[16]
18970
19306
19315

TV-RMOP model
19271
19271
19271

Table 4: The measures for test .

Initial
12135.5

Ratio
0.0455

RMSE
233.5

MAE
172.7

Table 5: Load in Shanghai region.

Year
1990
1991
1992
1993
1994
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010

Load (hundred million Kwh)

264.74
288.78
317.38
345.86
337.3
403.27
430.4
454.26
482.94
501.2
559.42
592.99
645.71
745.97
821.44
921.97
990.15
1072.38
1138.22
1153.38
1295.87

the fluctuation in fuzzy time series. Then the problem can be

formulated as
min

(initial, ratio)
= [RMSE, MAE, MAPE, MAP, THEIL]

subject to initial [12000, 13000] , ratio (0, 0.2] .

(5)

MAPE
0.9098

MAP
2.79

THEIL
0.01212

We use a MATLAB function gamultiobj to solve (5).

The initialization of genetic algorithm parameters is set as
follows. The fraction of population on non-dominated front
is 0.1, the number of individuals is 50, maximum number
of generations allowed is 100, and termination tolerance on
fitness function value is 1 100. The Pareto solutions of
training phase are provided in Table 2 (seen in Figure 3).
We select 12135.5 and 0.0455 as optimal parameters to
forecast the enrollment from 1990 to 1992. For this initial
and ratio, the universe of discourse is defined as [12135.5,
19694.1]. The intervals are finally defined as: 1 = [12135.5,
12681.6], 2 = [12681.6, 13252.3], 3 = [13252.3, 13848.6],
4 = [13848.6, 14471.8], 5 = [14471.8, 15123], 6 = [15123,
15803.6], 7 = [15803.6, 16514.7], 8 = [16514.7, 17257.9],
9 = [17257.9, 18034.5], 10 = [18034.5, 18846], 11 =
[18846, 19694.1]. Tables 3 and 4 list the forecast results and
corresponding measures, respectively. The best RMSE value
is 500 with the interval length of 500 for Chens method [20].
The RMSE value is 425.46 for Kuos model [15] and 327.54 for
Wongs model [16]. The proposed model produces 233.5 for
RMSE. Obviously, the TV-RMOP model gives more accurate
forecast than [15, 16, 20].

4. Applying TV-RMOP Model on Long Term

Load Forecasting
Long term load forecasting is needed for power supply and
delivery system. We use twenty years (19902010) load of
Shanghai region in China for forecast (data from the 2011
Shanghai statistical yearbook). The load is shown in Table 5.
Where loads from 1990 to 2007 are used for the training
(train ), and loads from 2008 to 2010 are used for test (test ).
The minimal value of train is 264.74, [240, 260] is taken
as constrains of variable initial. (0, 0.4] is taken as constrains

Mathematical Problems in Engineering

Table 6: Pareto solutions of load.

Initial
257.68
250.80
256.51
250.80
252.48

Ratio
0.2317
0.2400
0.2329
0.2400
0.2380

RMSE
45.28
44.06
45.22
44.06
44.24

MAE
34.47
34.49
34.54
34.49
34.30

MAPE
5.8185
5.9488
5.8476
5.9488
5.8942

MAP
12.184
13.627
12.448
13.627
13.273

THEIL
0.0719
0.0699
0.0718
0.0699
0.0702

Table 7: The forecasting results for test .

Year
2008
2009
2010

Actual load
1138.22
1153.38
1295.87

TV-RMOP model
1021.1
1266.2
1266.2

Table 8: The measures for test .

Initial
250.8

Ratio
0.24

RMSE
95.43

MAE
86.53

of the other variable ratio. Then the problem can be formulated as

min

MAPE
7.45

MAP
10.29

THEIL
0.07966

show that the TV-RMOP model is more accurate than existing models. Finally, we apply this model for long-term load
forecasting of Shanghai Region with satisfactory result.

(initial, ratio)
= [RMSE, MAE, MAPE, MAP, THEIL]

(6)

subject to initial [240, 260] , ratio (0, 0.4] .

The initialization of genetic algorithm parameters is set as
follows. The fraction of population on non-dominated front
is 0.1, the number of individuals is 50, maximum number
of generations allowed is 100, and termination tolerance on
fitness function value is 1 100. The Pareto solutions of
training phase are provided in Table 6, also seen in Figure 4.
Select 250.8 and 0.24 as optimal parameters to forecast the
loads from 2008 to 2010. Tables 7 and 8 list the forecast results
and corresponding measures, respectively. It is seen that the
MAPE and MAP measures are 7.45% and 10.29%, which meet
the actual demand.

5. Conclusions
Fuzzy time series forecasting method are getting quite popular in recent years. Most studies focus on the forecasting
rules and methods of defuzzification. How to determine the
lengths of intervals is less studied. The decision on what the
lengths will be is important for forecasting accuracy.
In this study, a new model TV-RMOP is proposed, which
formulates the forecasting as a multiobjective optimization
problem and solve it by the use of genetic algorithms. The
optimal lengths of intervals can be determined in the proposed model. The model is tested for forecasting the enrollment of University of Alabama. The experimental results

Acknowledgments
The authors would like to thank the anonymous reviewers in
MPE for helpful suggestions and correction. This work was
conducted by using of the MATLAB software and the authors
declared that it has no conflict of interests to this work.
The work was financially supported by Shanghai Municipal
Nature Science Foundation under Grant 10ZR1401400.

References
[1] E. Almeshaiei and H. Saltan, A methodology for electric power
load forecasting, Alexandria Engineering Journal, vol. 50, pp.
137144, 2011.
[2] H. Hahn, S. Meyer-Nieberg, and S. Pickl, Electric load forecasting methods: tools for decision making, European Journal
of Operational Research, vol. 199, no. 3, pp. 902907, 2009.
[3] Y. M. Wi, S. K. Joo, and K. B. Song, Holiday load forecasting
using fuzzy polynomial regression with weather feature selection and adjustment, IEEE Transactions on Power Systems, vol.
27, no. 2, pp. 596603, 2012.
[4] S. Fan and R. J. Hyndman, Short-term load forecasting based
on a semi-parametric additive model, IEEE Transactions on
Power Systems, vol. 27, no. 1, pp. 134141, 2012.
[5] J. Taylor, Short-term load forecasting with exponentially
weighted methods, IEEE Transactions on Power Systems, vol.
27, no. 1, pp. 458464, 2012.
[6] A. Khosravi, S. Nahavandi, D. Creighton, and D. Srinivasan,
Interval type-2 fuzzy logic systems for load forecasting: a comparative study, IEEE Transactions on Power Systems, vol. 27, no.
3, pp. 12741282, 2012.

Mathematical Problems in Engineering

[7] M. Xia, Y. C. Zhang, L. G. Wen, and X. L. Ye, Fashion retailing
forecasting based on extreme learning machine with adaptive
metrics of inputs, Knowledge-Based Systems, vol. 36, pp. 253
259, 2012.
[8] N. Amjady and A. Daraeepour, Midterm demand prediction of
electrical power systems using a new hybrid forecast technique,
IEEE Transactions on Power Systems, vol. 26, no. 2, pp. 755765,
2011.
[9] R. A. Hooshmand, H. Amooshahi, and M. Parastegari, A
hybrid intelligent algorithm based short-term load forecasting
approach, International Journal of Electrical Power & Energy
Systems, vol. 45, no. 1, pp. 313324, 2013.
[10] S. Zhaowei, Z. Lingfeng, M. Shangjun, F. Bin, and Z. Taiping,
Incomplete time series prediction using max-margin classification of data with absent features, Mathematical Problems in
Engineering, vol. 2010, Article ID 513810, 14 pages, 2010.
[11] E. Tian, D. Yue, and Z. Gu, Robust control for nonlinear
systems over network: a piecewise analysis method, Fuzzy Sets
and Systems, vol. 161, no. 21, pp. 27312745, 2010.
[12] E. G. Tian and D. Yue, Reliable filter design for T-S fuzzy
model-based networked control systems with random sensor
failure, International Journal of Robust and Nonlinear Control,
vol. 23, no. 1, pp. 1532, 2013.
[13] Q. Song and B. S. Chissom, Forecasting enrollments with fuzzy
time seriespart I, Fuzzy Sets and Systems, vol. 54, no. 1, pp. 1
9, 1993.
[14] Q. Song and B. S. Chissom, Forecasting enrollments with fuzzy
time seriespart II, Fuzzy Sets and Systems, vol. 62, no. 1, pp.
18, 1994.
[15] I. H. Kuo, S. J. Horng, T. W. Kao, T. L. Lin, C. L. Lee, and Y.
Pan, An improved method for forecasting enrollments based
on fuzzy time series and particle swarm optimization, Expert
Systems with Applications, vol. 36, no. 3, pp. 61086117, 2009.
[16] W. K. Wong, E. Bai, and A. W. C. Chu, Adaptive time-variant
models for fuzzy-time-series forecasting, IEEE Transactions on
Systems, Man, and Cybernetics Part B, vol. 40, no. 6, pp. 1531
1542, 2010.
[17] R. T. Marler and J. S. Arora, Survey of multi-objective optimization methods for engineering, Structural and Multidisciplinary Optimization, vol. 26, no. 6, pp. 369395, 2004.
[18] Y. Tang, Z. Wang, H. Gao, S. Swift, and J. Kurths, A constrained
evolutionary computation method for detecting controlling
regions of cortical networks, IEEE/ACM Transactions on Computational Biology and Bioinformatics, vol. 9, no. 6, pp. 1569
1581, 2012.
[19] Y. Tang, Z. Wang, W. K. Wong, J. Kurths, and J.-A. Fang, Multiobjective synchronization of coupled systems, Chaos, vol. 21,
no. 2, Article ID 025114, 12 pages, 2011.
[20] S. M. Chen, Forecasting enrollments based on fuzzy time
series, Fuzzy Sets and Systems, vol. 81, no. 3, pp. 311319, 1996.

Advances in

Operations Research
Hindawi Publishing Corporation
http://www.hindawi.com

Volume 2014

Advances in

Decision Sciences
Hindawi Publishing Corporation
http://www.hindawi.com

Volume 2014

Journal of

Applied Mathematics

Algebra

Hindawi Publishing Corporation

http://www.hindawi.com

Hindawi Publishing Corporation

http://www.hindawi.com

Volume 2014

Journal of

Probability and Statistics

Volume 2014

The Scientific
World Journal
Hindawi Publishing Corporation
http://www.hindawi.com

Hindawi Publishing Corporation

http://www.hindawi.com

Volume 2014

International Journal of

Differential Equations
Hindawi Publishing Corporation
http://www.hindawi.com

Volume 2014

Volume 2014

Submit your manuscripts at

http://www.hindawi.com
International Journal of

Advances in

Combinatorics
Hindawi Publishing Corporation
http://www.hindawi.com

Mathematical Physics
Hindawi Publishing Corporation
http://www.hindawi.com

Volume 2014

Journal of

Complex Analysis
Hindawi Publishing Corporation
http://www.hindawi.com

Volume 2014

International
Journal of
Mathematics and
Mathematical
Sciences

Mathematical Problems
in Engineering

Journal of

Mathematics
Hindawi Publishing Corporation
http://www.hindawi.com

Volume 2014

Hindawi Publishing Corporation

http://www.hindawi.com

Volume 2014

Volume 2014

Hindawi Publishing Corporation

http://www.hindawi.com

Volume 2014

Discrete Mathematics

Journal of

Volume 2014

Hindawi Publishing Corporation

http://www.hindawi.com

Discrete Dynamics in
Nature and Society

Journal of

Function Spaces
Hindawi Publishing Corporation
http://www.hindawi.com

Abstract and
Applied Analysis

Volume 2014

Hindawi Publishing Corporation

http://www.hindawi.com

Volume 2014

Hindawi Publishing Corporation

http://www.hindawi.com

Volume 2014

International Journal of

Journal of

Stochastic Analysis

Optimization

Hindawi Publishing Corporation

http://www.hindawi.com

Hindawi Publishing Corporation

http://www.hindawi.com

Volume 2014

Volume 2014