import pymks

%matplotlib inline
%load_ext autoreload
%autoreload 2

import numpy as np
import matplotlib.pyplot as plt


"""
Explanation: Effective Stiffness of Fiber Composite
Introduction
This example demonstrates the use of the homogenization model from pyMKS on a set of fiber-like structures.  These structures are simulated to emulate fiber-reinforced polymer samples.  For a summary of homogenization theory and its use with effective stiffness properties please see the Effective Siffness example.  This example will first generate a series of random microstructures with various fiber lengths and volume fraction.  The ability to vary the volume fraction is a new functionality of this example.  Then the generated stuctures will be used to calibrate and test the model based on simulated effective stress values.  Finally we will show that the simulated response compare favorably with those generated by the model.  
Generating Structures
These first lines inport important packages that will be used to run pymks.
End of explanation
"""


sample_size = 100
n_samples = 4 * [sample_size]
size = (101, 101)
elastic_modulus = (1.3, 75)
poissons_ratio = (0.42, .22)
macro_strain = 0.001
n_phases = 2
grain_size = [(40, 2), (10, 2), (2, 40), (2, 10)]
v_frac = [(0.7, 0.3), (0.6, 0.4), (0.3, 0.7), (0.4, 0.6)]
per_ch = 0.1


"""
Explanation: Now we are defining the parameters which we will use to create the microstructures.  n_samples will determine how many microstructures of a particular volume fraction we want to create.  size determines the number of pixels we want to be included in the microstructure.  We will define the material properties to be used in the finite element in elastic_modulus, poissons_ratio and macro_strain.  n_phases and grain_size will determine the physical characteristics of the microstructure.  We are using a high aspect ratio in creating our microstructures to simulate fiber-like structures.  The volume_fraction variable will be used to vary the fraction of each phase.  The sum of the volume fractions must be equal to 1.  The percent_variance variable introduces some variation in the volume fraction up to the specified percentage.
End of explanation
"""


from pymks.datasets import make_elastic_stress_random


dataset, stresses = make_elastic_stress_random(n_samples=n_samples, size=size, grain_size=grain_size,
                                                 elastic_modulus=elastic_modulus, poissons_ratio=poissons_ratio,
                                                 macro_strain=macro_strain, volume_fraction=v_frac,
                                                 percent_variance=per_ch)


"""
Explanation: Now we will create the microstructures and generate their responses using the make_elastic_stress_random function from pyMKS.  Four datasets are created to create the four different volume fractions that we are simulating.  Then the datasets are combined into one variable.   The volume fractions are listed in the variable v_frac.  Variation around the specified volume fraction can be obtained by varying per_ch.  The variation is randomly generated according a uniform distribution around the specified volume fraction.
End of explanation
"""


from pymks.tools import draw_microstructures
examples = dataset[::sample_size]
draw_microstructures(examples)


"""
Explanation: Now we are going to print out a few microstructres to look at how the fiber length, orientation and volume fraction are varied.
End of explanation
"""


from pymks import MKSHomogenizationModel
from pymks import PrimitiveBasis


p_basis = PrimitiveBasis(n_states=2, domain=[0, 1])
model = MKSHomogenizationModel(basis=p_basis, correlations=[(0, 0)], periodic_axes=[0, 1])


"""
Explanation: Creating the Model
Next we are going to initiate the model. The MKSHomogenizationModel takes in microstructures and runs two-point statistics on them to get a statistical representation of the microstructures.  An expalnation of the use of two-point statistics can be found in the Checkerboard Microstructure Example.  Then the model uses PCA and regression models to create a linkage between the calcualted properties and structures.
Here we simply initiate the model.
End of explanation
"""


from sklearn.cross_validation import train_test_split


flat_shape = (dataset.shape[0],) + (dataset[0].size,)
data_train, data_test, stress_train, stress_test = train_test_split(
    dataset.reshape(flat_shape), stresses, test_size=0.2, random_state=3)


"""
Explanation: Now we are going to split our data into testing and training segments so we can test and see if our model can accurately predict the effective stress.
End of explanation
"""


from sklearn.grid_search import GridSearchCV


params_to_tune = {'degree': np.arange(1, 4), 'n_components': np.arange(2, 8)}
fit_params = {'size': dataset[0].shape}
gs = GridSearchCV(model, params_to_tune, fit_params=fit_params).fit(data_train, stress_train)


"""
Explanation: We will use sklearn's GridSearchCV to optimize the n_components and degree for our model. Let's search over the range of 1st order to 3rd order polynomial for degree and 2 to 7 principal components for n_components.
End of explanation
"""


print('Order of Polynomial'), (gs.best_estimator_.degree)
print('Number of Components'), (gs.best_estimator_.n_components)
print('R-squared Value'), (gs.score(data_test, stress_test))


from pymks.tools import draw_gridscores_matrix

draw_gridscores_matrix(gs, ['n_components', 'degree'], score_label='R-Squared',
                       param_labels=['Number of Components', 'Order of Polynomial'])


from pymks.tools import draw_gridscores

gs_deg_1 = [x for x in gs.grid_scores_ \
            if x.parameters['degree'] == 1]
gs_deg_2 = [x for x in gs.grid_scores_ \
            if x.parameters['degree'] == 2]
gs_deg_3 = [x for x in gs.grid_scores_ \
            if x.parameters['degree'] == 3]

draw_gridscores([gs_deg_1,  gs_deg_2, gs_deg_3], 'n_components', 
                data_labels=['1st Order', '2nd Order', '3rd Order'],
                param_label='Number of Components', score_label='R-Squared')


"""
Explanation: Let's take a look at the results.
End of explanation
"""


model = gs.best_estimator_

"""
Explanation: Our best model was found to have degree equal to 3 and n_components equal to 5. Let's go ahead and use it.
End of explanation
"""


from pymks.tools import draw_components_scatter


stress_predict = model.predict(data_test)
draw_components_scatter([model.reduced_fit_data[:, :3],
                         model.reduced_predict_data[:, :3]],
                        ['Training Data', 'Testing Data'],
                        legend_outside=True)


"""
Explanation: Structures in PCA space
Now we want to draw how the samples are spread out in PCA space and look at how the testing and training data line up.
End of explanation
"""


from pymks.tools import draw_goodness_of_fit


fit_data = np.array([stresses, model.predict(dataset)])
pred_data = np.array([stress_test, stress_predict])
draw_goodness_of_fit(fit_data, pred_data, ['Training Data', 'Testing Data'])


"""
Explanation: It looks like there is pretty good agreement between the testing and the training data.  We can also see that the four different fiber sizes are seperated in the PC space.  
Draw Goodness of fit
Now we are going to look at how well our model predicts the properties of the structures.  The calculated properties will be plotted against the properties generated by the model.  We should see a linear realtionship with a slope of 1.
End of explanation
"""


