//======    
int n;
//======  
double U0, UN;

//======   
//======  
double f()
{
	//======   
	static int raw = -1, k = -1, col = 0;

	//======   
	col++;

	//====== k   
	//======    
	if (++k % n == 0)
	{
		col = 0;		//   
		raw++;			//   
	}

	//======   
	return col==raw ? -2. 
		  : col == raw-1 || col==raw+1 ? 1. 
		  : 0.;
}

double fu()
{
	//====       (5)
	static double
		dU = (UN-U0)/(n+1),
		d = U0;
	return d += dU;
}



void main()
{
	//=======     
	n = 4;
	U0 = 100.;
	UN = 0.;

	//=======  valarray ( )
	int nn = n*n;

	//=======    
	valarray<double> a(nn), u(n), v(n);
	
	//=======   
	generate (&a[0], &a[nn], f);
	generate (&u[0], &u[n], fu);

	out("Initial matrix", a);
	out("Initial vector", u);
	
	//=======    
	for (int i=0; i<n; i++)
	{
		//=======  i-  
		valarray<double> s = a[slice (i*n, n , 1)];

		//=======    
		//=======     v
		transform(&s[0], &s[n], &u[0], &v[0],
						multiplies<double>());

		//=======  , 
		//======= i-    
		cout << "\nb[" << i << "] = " << v.sum();
	}

	cout<<"\n\n";
}


for (int i=0, id=0;  i<n;  i++, id+=n)
{
	transform(&a[id], &a[id+n], &u[0], &v[0],
											multiplies<double>());
	cout << "\nb[" << i << "] = " << v.sum();
}



void CChildView::OnPaint() 
{
	CPaintDC dc(this);
	CGraph(m_Points, "Field Distribution",
							"x[m]","Field").Draw(&dc);

}



#pragma once
#include "Graph.h"

class CChildView : public CWnd
{
	//=====     
	friend class CParamDlg;
	friend class CGraph;
private:
	//=====    
	vector<CDPoint> m_Points;
	//=====      (. f  ()
	vector<double> m_f, m_r;

	//=====   (. N)
	int m_n;
	//===== 
	double	m_k,			//  k
			m_L,			//   
			m_g0,			// ,   
			m_d0,
			m_gn,			// ,   
			m_dn;
	CParamDlg *m_pDlg;	//   

public:
	CChildView();
	virtual ~CChildView();

	virtual BOOL PreCreateWindow(CREATESTRUCT& cs);

	//=====   
	void Resize();
	//=====    
	void Solve();
protected:
	afx_msg void OnPaint();
	DECLARE_MESSAGE_MAP()
};



CChildView::CChildView()
{
	m_n = 200;
	m_k = -0.0005;
	m_L = 200.;
	//======     Uo=100
	m_g0 = 0.;
	m_d0 = 100.;
	m_gn = 0.;
	m_dn = 0.;
	
	Resize();
	m_pDlg = 0;
}



CChildView::~CChildView()
{
	m_Points.clear();
	m_f.clear();
	m_r.clear();
}



void CChildView::Resize()
{
	//=====    N+1 (  0- )
	int n = m_n + 1;

	m_Points.resize(n, CDPoint(0.,0.));
	m_f.resize(n, 0.);
	m_r.resize(n, 1.);
} 



void CChildView::Solve()
{	
	Resize();
	int n = m_n + 1;

	//=======   
	vector<double> a(n), b(n), c(n);
	
	//=======  
	vector<double> d(n), e(n);

	double	h = m_L / m_n,		//     x
			hh = h * h;			//  

	//=======   0-   
	a[0] = 0.;
	b[0] = 0.;
	c[0] = 0.;

	//=======   x   
	m_Points[0].x = 0.;

	for (int i=1;  i < m_n;  i++)
	{
		m_Points[i].x = i * h;
	//=======   (4)
		a[i] = m_r[i-1]/hh;
		c[i] = m_r[i]/hh;
		b[i] = - a[i] - c[i] + m_k;
	}
	m_Points[m_n].x = m_L;

	//=======   
	d[0] = m_g0;		//  
	e[0] = m_d0;
	double den;

	for (i=1; i < m_n; i++)
	{
		//=======  
		den = a[i] * d[i-1] + b[i];
		d[i] = -c[i] / den;
		e[i] = (m_f[i] - a[i] * e[i-1]) / den;
	}

	//=======    
	den = 1. - m_gn * d[m_n-1];

	//=======    
	if (den==0.)
	{
		MessageBox("  ","",MB_OK);
		return;
	}

	//=======      
	//=======   (13)
	m_Points[m_n-1].y = (e[m_n-1] + m_dn * d[m_n-1])/den;
	m_Points[m_n].y = (m_dn + m_gn* e[m_n-1])/den;

	//=======   
	for (i = m_n-2;  i >= 0;  i--)
		m_Points[i].y = d[i] * m_Points[i+1].y + e[i];

	Invalidate();
}



int CChildView::OnCreate(LPCREATESTRUCT lpCreateStruct)
{
	if (CWnd::OnCreate(lpCreateStruct) == -1)
		return -1;
	//=======  ,   
	Solve();
	return 0;
}



#pragma once

class CDPoint
{
public:
	//=======      
	double x, y;

	//=======     
	CDPoint()
	{
		x=0.;  y=0.;
	}
	CDPoint(double xx, double yy)
	{
		x=xx;  y=yy;
	}

	CDPoint& operator=(const CDPoint& pt)
	{
		x = pt.x;
		y = pt.y;
		return *this;
	}

	CDPoint(const CDPoint& pt)
	{
		*this = pt;
	}
};

	//=====  , 
	//=====      
struct TData
{
	//=====     
	int Power;
	//=====   X
	bool bX;
	double
		//======= 
		Min, Max,
		//=======  (10   Power)
		Factor,
		//=======    ()
		Step,
		//=======  
		dStep,
		//=======     ()
		Start, End,
		//=======    
		dStart, dEnd;
};

	//===== ,    
class CGraph
{
public:
	//===== ,    
	TData m_DataX, m_DataY;
	//=====   
	vector <CDPoint>& m_Points;
	//=====    
	CSize m_Size;
	//=====    
	CPoint m_Center;
	//=====    
	CString m_sTitle, m_sX, m_sY;
	//=====   
	CPen m_Pen;
	//=====   
	CFont m_TitleFont, m_Font;
	//=====   (  )
	int		m_LH,
	//=====  
			m_Width;
	//=====  
	COLORREF m_Clr;

	//=======    
	CGraph(vector<CDPoint>& pt, CString sTitle,
					CString sX, CString sY);
	virtual ~CGraph();
	//=====  TData    
	void Scale(TData& data);
	//=====     
	int MapToLogX (double d);
	int MapToLogY (double d);
	//=====    
	void Draw (CDC *pDC);
	//=====   
	void DrawLine(CDC *pDC);
	//=====     
	CString MakeLabel(bool bX, double& d);
};



#include "StdAfx.h"
#include "graph.h"

//=====  ,  
#define SCALE_X 0.6
#define SCALE_Y 0.6

//=====      
void gScale (double span, double& step)
{
	//=====  span   
	//=====      
	//=====   ,  
	int power = int(floor(log10(span)));
	//=====  (zoom factor)
	double factor = pow(10, power);
	//=====   ( 1 < span < 10)
	span /= factor;

	//=====    
	if (span<1.99)
		step=.2;
	else if (span<2.49)
		step=.25;
	else if (span<4.99)
		step=.5;
	else if (span<10.)
		step= 1.;

	//=====     (step*10^power)
	step *= factor; 
}



void CGraph::Scale (TData& data)
{
	//=====     
	if (m_Points.empty())
		return;

	//=====   
	data.Max = data.bX ? m_Points[0].x : m_Points[0].y;
	data.Min = data.Max;

	//=====  	
	for (UINT j=0;  j<m_Points.size();  j++)
	{
		double d = data.bX ? m_Points[j].x 
						   : m_Points[j].y; 
		if (d < data.Min)
			data.Min = d;
		if (d > data.Max)
			data.Max = d;
	}

	//=====     
	double ext = max(fabs(data.Min),fabs(data.Max));

	//=====    
	//=====  3 ,      7 ,
	//=====       
	double power = ext > 0.? log10(ext) + 3. : 0.;	
	data.Power = int(floor(power/7.));

	//=====       
	if (data.Power != 0)
		//=====     
		data.Power = int(floor(power)) - 3;
	//=====  
	data.Factor = pow(10,data.Power);	

	//=====   
	double span = (data.Max - data.Min)/data.Factor;
	//=====   ,
	if (span == 0.)
		span = 0.5; //    

	//=====      
	gScale (span, data.Step);

	//=====     
	data.dStep = data.Step * data.Factor;

	//=====       
	//=====    
	data.dStart = data.dStep * int(floor(data.Min/data.dStep));
	data.Start = data.dStart/data.Factor;

	//=====    
	for (data.End = data.Start;
		data.End < data.Min/data.Factor + span-1e-10;  
		data.End += data.Step)
		;
	data.dEnd = data.End*data.Factor;
}



//=======   CGraph
CGraph::CGraph (vector<CDPoint>& pt, CString sTitle,
					CString sX, CString sY)
					: m_Points(pt)
{
	//=======  ,   	ZeroMemory(&m_DataX, sizeof(TData));
	ZeroMemory(&m_DataY, sizeof(TData));
	m_DataX.bX = true;
	m_DataY.bX = false;
	m_sTitle = sTitle;
	m_sX = sX;
	m_sY = sY;

	//=======     
	m_Font.CreateFont(16,0,0,0,100,0,0,0,DEFAULT_CHARSET,
				OUT_RASTER_PRECIS,CLIP_DEFAULT_PRECIS,
				DEFAULT_QUALITY,FF_DONTCARE,"Arial");
	//=======    
	TEXTMETRIC tm;
	CClientDC dc(0);
	dc.SelectObject(&m_Font);
	dc.GetTextMetrics(&tm); 
	m_LH = tm.tmHeight;

	//=======     
	m_TitleFont.CreateFont(24,0,0,0,100,0,0,0,DEFAULT_CHARSET,
				OUT_RASTER_PRECIS, CLIP_DEFAULT_PRECIS,
				DEFAULT_QUALITY,FF_DONTCARE,"Times New Roman");

	//=======   
	m_Clr = RGB(0,0,255);
	m_Width = 2;
}



int CGraph::MapToLogX (double d)
{
	return m_Center.x + int (SCALE_X * m_Size.cx * d);
}

int CGraph::MapToLogY (double d)
{
	return m_Center.y - int (SCALE_Y * m_Size.cy * d);
}

//=======  
CGraph::~CGraph(){}



void CGraph::Draw(CDC *pDC)
{
	//======    
	//======   ,  
	CWnd *pWnd = pDC->GetWindow();
	CRect r;
	pWnd->GetClientRect(&r);

	//======   
	m_Size = r.Size();
	m_Center = CPoint(m_Size.cx/2, m_Size.cy/2);
	
	int nDC = pDC->SaveDC();
	
	//======      
	CPen pen(PS_SOLID, 0, COLORREF(0));
	pDC->SelectObject(&pen);

	//======   
	int lt = MapToLogX(-0.5),
		rt = MapToLogX(0.5),
		tp = MapToLogY(0.5),
		bm = MapToLogY(-0.5);

	pDC->Rectangle (lt, tp, rt, bm);

	//======     
	pDC->SetTextColor(0);
	pDC->SetTextAlign(TA_LEFT | TA_BASELINE);
	
	//======  
	pDC->SelectObject (&m_Font);
	
	//======    
	Scale(m_DataX);
	Scale(m_DataY);

	//======   	
	CString s;
	s.Format("Min = %.3g",m_DataY.Min);
	pDC->TextOut(rt+m_LH, tp+m_LH, s);

	s.Format("Max = %.3g",m_DataY.Max);
	pDC->TextOut(rt+m_LH, tp+m_LH+m_LH, s);
		
	//======    
	CPen gridPen(PS_SOLID, 0, RGB(92,200,178));
	pDC->SelectObject(&gridPen);
	pDC->SetTextAlign(TA_CENTER | TA_BASELINE);
	
	//======    
	for (double x = m_DataX.Start;  
			x < m_DataX.End - m_DataX.Step/2.;
			x += m_DataX.Step)
	{
		//======   x
		double	xn = (x - m_DataX.Start) /
						(m_DataX.End - m_DataX.Start) - 0.5;

		//======   
		int xi = MapToLogX(xn);	
		//======   ,
		//======      
		if (x > m_DataX.Start && x < m_DataX.End)
		{
			pDC->MoveTo(xi, bm);
			pDC->LineTo(xi, tp);
		}
		//======   
		pDC->TextOut(xi, bm+m_LH, MakeLabel(true, x));
	}

	//=====      
	pDC->SetTextAlign(TA_RIGHT | TA_BASELINE);
	for (double y = m_DataY.Start;
			y < m_DataY.End - m_DataY.Step/2.;
			y += m_DataY.Step)
	{
		double yn = (y - m_DataY.Start) /
						(m_DataY.End - m_DataY.Start) - 0.5;

		int yi = MapToLogY(yn);
		if (y > m_DataY.Start && y < m_DataY.End)
		{
			pDC->MoveTo(lt, yi);
			pDC->LineTo(rt, yi);
			pDC->TextOut(lt-m_LH/2, yi, MakeLabel(false, y));
		}
	}

	//======   
	pDC->TextOut(lt-m_LH/2, tp - m_LH, m_sY);
	pDC->SetTextAlign(TA_LEFT | TA_BASELINE);
	pDC->TextOut(rt-m_LH, bm + m_LH, m_sX);

	//======  
	if (m_sTitle.GetLength() > 40)
		m_sTitle.Left(40);
	pDC->SelectObject(&m_TitleFont);
	pDC->SetTextAlign(TA_CENTER | TA_BASELINE);
	pDC->TextOut((lt+rt)/2, tp - m_LH, m_sTitle);
	
	//======   
	DrawLine(pDC);
	//======   GDI
	pDC->RestoreDC(nDC);
}



void CGraph::DrawLine(CDC *pDC)
{
	//======   
	if (m_Pen.m_hObject)
		m_Pen.DeleteObject();
	//======  
	m_Pen.CreatePen(PS_SOLID, m_Width, m_Clr);

	pDC->SelectObject(&m_Pen);

	double	x0 = m_DataX.dStart,
			y0 = m_DataY.dStart,
			dx = m_DataX.dEnd - x0,
			dy = m_DataY.dEnd - y0;

	for (UINT i=0;  i < m_Points.size();  i++)
	{
		//======  
		double	x = (m_Points[i].x - x0) / dx - .5,
				y = (m_Points[i].y - y0) / dy - .5;

		//======    
		CPoint pt (MapToLogX(x),MapToLogY(y));
		//======   ,   
		if (i==0)
			pDC->MoveTo(pt);
		else
			pDC->LineTo(pt);
	}
}



CString CGraph::MakeLabel(bool bX, double& v)
{
	CString s = "0.0";
	if (v == 0.)
		return s;

	//======    
	//======    20 
	s.Format("%20.10f",v);
	//======    ,
	//======    ()
	int nDigits = int(ceil(-log10(bX ? m_DataX.Step
									: m_DataY.Step)));
	//======      ,
	//======       
	if (nDigits <= 0)
		nDigits = -1;
	else
		if (bX)
			nDigits++;	//  

	//======   
	s.TrimLeft();

	//======    
	s = s.Left(s.Find(".") + nDigits + 1);

	int iPower = bX ? m_DataX.Power : m_DataY.Power;
	//======   ?
	if (iPower != 0)
	{
		//====== ,     (10^-3, 10^+4)
		CString add;
		add.Format("e%+d",iPower);
		s += add;
	}
	return s;
}




IDD_PARAM
  Source:
IDC_SOURCE
  Start  Source:
IDC_SOURCE1
  End  Source:
IDC_SOURCE2
  Value:
IDC_PROP
  Start  Properties:
IDC_PROP1
  End  Properties:
IDC_PROP2
  Nodes:
IDC_NODES
  Distance:
IDC_DIST
  Decrement:
IDC_DECR
  g  Left Boundary:
IDC_LEFTG
  d  Left Boundary:
IDC_LEFTD
  g  Right Boundary:
IDC_RIGHTG
  d  Right Boundary:
IDC_RIGHTD
 Add  Source:
IDC_ADDSOURCE
 Add  Properties:
IDC_ADDPROP
 Apply:
IDC_APPLY
 Close
IDCANCEL



#pragma once

class CParamDlg : public CDialog
{
	//=====    
	friend class CChildView;
	DECLARE_DYNAMIC(CParamDlg)

public:
	//=====    
	CChildView *m_pView;
	//=====     
	CParamDlg(CChildView* p);
	virtual ~CParamDlg();

	// Dialog Data
	enum { IDD = IDD_PARAM };

protected:
	virtual void DoDataExchange(CDataExchange* pDX);

	DECLARE_MESSAGE_MAP()
}; 



//====   
double	m_Source;
//====   ,    
int		m_SrcId1;
//====   ,    
int		m_SrcId2;
//====     
double	m_Prop;
//====        m_Prop
int		m_PropId1;
int		m_PropId2;



void CParamDlg::DoDataExchange(CDataExchange* pDX)
{
	DDX_Text(pDX, IDC_PROP2, m_PropId2);
	DDX_Text(pDX, IDC_PROP1, m_PropId1);
	DDX_Text(pDX, IDC_PROP, m_Prop);
	DDX_Text(pDX, IDC_SOURCE2, m_SrcId2);
	DDX_Text(pDX, IDC_SOURCE1, m_SrcId1);
	DDX_Text(pDX, IDC_SOURCE, m_Source);
	
	//=========     
	DDX_Text(pDX, IDC_NODES, m_pView->m_n);
	DDX_Text(pDX, IDC_DIST, m_pView->m_L);
	DDX_Text(pDX, IDC_DECR, m_pView->m_k);
	DDX_Text(pDX, IDC_LEFTG, m_pView->m_g0);
	DDX_Text(pDX, IDC_LEFTD, m_pView->m_d0);
	DDX_Text(pDX, IDC_RIGHTG, m_pView->m_gn);
	DDX_Text(pDX, IDC_RIGHTD, m_pView->m_dn);

	CDialog::DoDataExchange(pDX);
}



void CParamDlg::OnClickedApply(void)
{
	//======    
	UpdateData();
	//======      
	m_pView->Solve();
}

void CParamDlg::OnClickedAddsource(void)
{
	UpdateData();
	//======   m_f ( )
	for (int i=m_SrcId1;  i <= m_SrcId2;  i++)
	{
		if (0 <= i && i < m_pView->m_n)
			m_pView->m_f[i] = -m_Source;
	}
	m_pView->Solve();
}

void CParamDlg::OnClickedAddprop(void) 
{
	UpdateData();
	//======   m_r ( )
	for (int i=m_PropId1;  i <= m_PropId2;  i++)
	{
		if (0 <= i && i < m_pView->m_n && m_Prop > 0.)
			m_pView->m_r[i] = m_Prop;
	}
	m_pView->Solve();
}

void CParamDlg::OnClickedCancel(void)
{
	//======   
	m_pView->m_pDlg = 0;
	DestroyWindow();
}



#include "stdafx.h"
#include "Heat.h"
#include "ParamDlg.h"

IMPLEMENT_DYNAMIC(CParamDlg, CDialog)

CParamDlg::CParamDlg(CChildView* p)
	: CDialog(CParamDlg::IDD, p)
{
	m_pView = p;
	//=====    
	//=====    
	m_Prop = 1.0;
	m_PropId1 = 0;
	m_PropId2 = 0;
	m_Source = 0.0;
	m_SrcId1 = 0;
	m_SrcId2 = 0;
}

CParamDlg::~CParamDlg()
{
}



BOOL CParamDlg::OnInitDialog(void)
{
	CDialog::OnInitDialog();
	CRect r;
	//=====    
	//=====    
	CClientDC dc(this);
	int w = dc.GetDeviceCaps(HORZRES);
	int h = dc.GetDeviceCaps(VERTRES);

	//=====    
	GetWindowRect(&r);
	//=====     
	r.OffsetRect(w-r.right-10,h-r.bottom-30);
	MoveWindow(&r);
	return TRUE;
}



void CChildView::OnEditParameters(void)
{
	//=====    ,
	if (!m_pDlg)
	{
		//=====     
		m_pDlg = new CParamDlg(this);
		//=====       
		m_pDlg->Create(IDD_PARAM);
	}
}



void CParamDlg::PostNcDestroy(void)
{
	delete this;
}

