北航朝花夕拾bt
篇一：北航数值分析大作业第八题
北 京 航 空 航 天 大 学
数值分析大作业八
学院名称 自动化 专业方向 控制工程 学 号 学生姓名 许阳 教
师 孙玉泉
日 期 2014 年 11月26 日
一．题目
关于x , y , t , u , v , w 的方程组(A.3)
?0.5cost?u?v?w?x?2.67?t?0.5sinu?cosv?w?y?1.07?
(A.3) ?
0.5t?u?cosv?w?x?3.74???t?0.5u?v?sinw?y?0.79
以及关于z , t , u 的二维数表（见表A-1）确定了一个二元函数z=f(x , y)。
表A-1 二维数表
（ x,y)|0?x?0.8,0.5?y?1.5}上的1. 试用数值方法求出f (x , y) 在区域D?{
近似表达式
p(x,y)???crsxrys
i?0j?0
kk
要求p(x , y)以最小的k值达到以下的精度
????[f(xi,yi)?p(xi,yi)]2?10?7
i?0j?0
1020
其中xi?0.08i,yi?0.5?0.05j。
**
2. 计算f(xi*,y*j),p(xi,yj) (i=1,2,…,8 ; j=1,2,…,5) 的值，以观察p(x , y) 逼
近f (x , y)的效果，其中xi*?0.1i,y*j?0.5?0.2j。
二．算法设计
（一）总体思路
1.题目要求p(x,y)???crsxrys对f(x, y) 进行拟合，可选用乘积型最小二乘
i?0j?0k
k
拟合。(xi,yi)与f(xi,yi)的数表由方程组与表A-1得到。
**2.f(xi*,y*j)与1使用相同方法求得，p(xi,yj)由计算得出的p(x,y)直接带入
(xi*,y*j)求得。
（二）算法实现
1. (xi,yi)与f(xi,yi)的数表的获得
（ x,y)|0?x?0.8,0.5?y?1.5}上的f (x , y)值可由方程组及二维对区域D?{
数表得到。将区域D上的(xi,yi)分别回代于方程组(A.3)，成为关于t,u,v,w的4元非线性方程组，解出每个(xi,yj)对应的t,u。再通过表A-1进行插值近似，得到相应的z值。对应的z即为D区域上(xi,yj)对应的f(xi,yj)。从而得到(xi,yj)与f(xi,yi)的数表。
(1) 4元非线性方程组求解
(xi,yi)代入(A.3)后，原方程组变为关于t,u,v,w的4元非线性方程组。观察到方程组中方程形式较为简单，易于对变量t,u,v,w求偏导数，故而选用Newton法对方程组求解。
计算方程组矩阵为：
?0.5cost?u?v?w?xi?2.67?
???t?0.5sinu?cosv?w?yi?1.07?
F(t,u,v,w)??
0.5t?u?cosv?w?xi?3.74????t?0.5u?v?sinw?y?0.79?
i??
计算方程组偏导数矩阵为：
111???0.5sint
??
10.5cosu11??
F'(t,u,v,w)??
0.51?sinv1????10.51cosw???
迭代公式为：
?t(k?1)??t(k)???t(k)?
?(k?1)??(k)??(k)??u??u???u?
?v(k?1)???v(k)????v(k)?，k=0,1,2,…,n ???????w(k?1)??w(k)???w(k)???????
其中
??t(k)???t(k)?
?(k)??(k)?
???u?(k)(k)(k)(k)??u(k)(k)(k)(k)
F'(t,u,v,w)??F(t,u,v,w)的解。 为线性方程组??v(k)??(k)?
??v?????w(k)???w(k)?????
取max(|?t(k)|,|?u(k)|,|?v(k)|,|?w(k)|)/max(|t(k)|,|u(k)|,|v(k)|,|w(k)|)??为迭代终止条件。
由表A-1观察到t,u基本在（0,2）上，于是选取（t(0),u(0),v(0),w(0))?(1,1,1,1)为迭代初值。
通过以上方法求得与(xi,yj)对应的(tij,uij)。 (2) 分片二元双二次代数插值
为保证代数插值的收敛性，应采用分片低次插值。故此使用分片双二次代数插值。
tm?0.2m,un?0.4n(m?0,1,...,5;n?0,1,...,5)
给定(tij,uij)如满足如下关系式：
tm?0.1?tij?tm?0.1，1?m?4
un?0.2?uij?un?0.2，1?n?4
则选择(tk,ur)(k?m,m?1,m?2;r?n,n?1,n?2)为插值节点，相应插值多项式为
p22(tij,uij)???lk(tij)lr(uij)z(tk,ur)
k?mr?nm?2n?2
其中
lk(tij)??
m?2
tij?tqt?tq
（k?m,m?1，m?2)
q?mk
q?kn?2
lr(uij)??
q?nq?k
uij?uquk?uq
（r?n,n?1,n?2)
如果tij?0.3或tij?0.7，则上式取m=1或m=4；如果uij?0.6或uij?1.4，则取n=1或n=4。
得到p22(tij,uij)表达式后，直接带入(tij,uij)，得到的值即为与(xi,yj)对应的
f(xi,yj)。
2. 乘积型最小二乘曲面拟合
使用乘积型最小二乘拟合，根据k值不用，有基函数矩阵如下：
0k0k?x0??y0??x0?y0????B?????? ， G??????
?x0?xk??y0?yk?iij????j
数表矩阵如下：
?f(x0,y0)??U????
?f(x,y)?
i0?
f(x0,yj)?
??? f(xi,yj)??
记C=[crs]，则系数crs的表达式矩阵为：
-1T
C?（BTB）BUG(GTG)?1
通过求解如下线性方程，即可得到系数矩阵C。
(BTB)C(GTG)?BTUG
篇二：北航北海大学语文期末重点
大学语文重点
1、 先秦文学形态：诗歌 散文
2、 《诗经》是我国第一部诗歌总集，是我国诗歌传统的源头。
3、 历史散文代表作：《尚书》、《春秋》、《左传》、《国语》、《战国策》。
4、 诸子散文代表作：儒家《论语》、《孟子》。《荀子》，道家《老子》、《庄子》，法家《管子》、
《商君书》、《韩非子》。
《论语》是记录孔子及其弟子的言论和活动的语录体散文著作。由孔子弟子及再传弟子编
纂而成。
5、《楚辞》起源于战国时代的楚国，它是由我国伟大的浪漫主义爱国诗人屈原在楚国民歌基础
上创作的一种新体诗歌。 代表作《离骚》。
6、秦代文学 吕不韦《吕氏春秋》 全书分为十二记、八览、六论，称为“吕览”。
李斯《谏逐客书》
7、汉赋 枚乘的《七发》 它的出现为汉大赋的形成奠定了基础。
8、《史记》是司马迁创作的我国第一部纪传体通史，记载了中国历史上自黄帝、下至西汉武帝
时代3000年社会、历史和制度的兴衰沿革。
9、汉乐府是继《诗经》、《楚辞》之后，中国古代诗歌史上的一种新体诗。
10、古诗十九首是汉代文人五言诗中最成熟的作品。
11、建安文学是指以“三曹”和“七子”为主体的文学创作。
12、鲍照、谢灵运和颜延之并称为“元嘉三大家”。
13、庾信《拟咏怀》、《哀江南赋》
庾信的成就使他成为集汉魏南北朝诗歌艺术之大成的一代大家。
14、《西洲曲》是南朝乐府民歌中的杰作。
15、《木兰诗》是北朝乐府民歌中最杰出的作品。它与《孔雀东南飞》并称为我国古代诗歌史上
的“乐府双壁”。
16、志怪小说 干宝的《搜神记》 志人小说 刘义庆的《世说新语》
17、“初唐四杰”包括王勃、杨炯、卢照邻、骆宾王。
18、山水田园诗派 王维《夜归鹿门歌》 孟浩然《使至塞上》
19、边塞诗派 高适《燕歌行》 岑参《白雪歌送武判官归京》
王昌龄被誉为“七绝圣手”《从军行》
20、李白诗歌的艺术特色
（1）塑造了一个以才气自负，不同凡响，具有炽热情感、强烈个性的抒情主人公形象；题材选
择神异奇特、可敬可愕。
（2）奇思遐想与自然天真相结合。
（3）想象与比兴、象征和拟人、夸张有机地结合。
21、杜甫诗歌的艺术特色：沉郁顿挫。
22、古文运动：指从贞元到元和的二三十年间，由韩愈首倡把奇句单行、上承先秦两汉文体的
散文称为古文，并使之和六朝以来流行的骈文对立；后又得柳宗元的大力支持，古文的业绩更
著，影响更大，结果古文逐渐压倒了骈文，成为文坛的主要风尚。
韩愈《上襄阳于相公书》 宗元《永州八记》
23、以爱情为主题的作品如《任氏传》、《柳毅传》、《霍小玉传》、《李娃传》、《莺莺传》等，在
唐传奇中成就最高。
24、词是一种配合音乐歌唱的新型格律诗体。
25、温庭筠 花间派词人 《菩萨蛮》 李煜 南唐 《浣溪沙》
26、杨万里、范成大、陆游和尤袤号称“中兴四大诗人”。
27、话本原是说书艺人的底本，是随着民间“说话”伎艺发展起来的一种文学形式。
28、元散曲的主要作家
关汉卿 《窦娥冤》 马致远“秋思之祖”的《天净沙》
睢景臣《高祖还乡》 张养浩《潼关怀古》
29、四大名著
罗贯中《三国演义》 施耐庵《水浒传》 吴承恩《西游记》 曹雪芹《红楼梦》
30、汤显祖 《牡丹亭》“临川四梦”：《紫钗记》 《邯郸记》 《南柯记》《牡丹亭》
31、清初爱国诗人：顾炎武 黄宗羲 王夫之
32、清初戏曲：李玉 四种曲《一捧雪》 《人兽关》 《永团圆》 《占花魁》
李渔 《闲情偶寄》
33、洪昇 《长生殿》 孔尚任 《桃花扇》
34、乾嘉诗派诗人： 郑德潜 郑燮《家书》 袁枚《随园诗话》
35、通城派古文：是清中叶最著名的一个散文流派。主要作家有方苞、刘大櫆、姚鼐，他们都
是安徽桐城人，“桐城派”即因此得名。
方苞《狱中杂记》 刘大櫆《论文偶记》 姚鼐《登泰山记》
36、鲁迅《呐喊》、《彷徨》
《狂人日记》是中国现代文学史上的第一篇小说。 郁达夫的《沉沦》
37、诗歌：胡适 他的《尝试集》是最早出版的一部新诗集。
刘半农的《教我如何不想她》 郭沫若：《女神》 《凤凰涅磐》
38、“新月诗派”即“新格律诗派”，以闻一多和徐志摩为代表人物的自由诗派。
闻一多《红烛》和《死水》 徐志摩《志摩的诗》
“新月诗派”核心主张是讲究诗的“三美”：音乐美、绘画美、建筑美。
39、鲁迅的散文诗《野草》 散文集《朝花夕拾》
40、小说 茅盾《子夜》 巴金《激流三部曲》 尤其是《家》
老舍《四世同堂》 《骆驼祥子》
41、京派小说 沈从文的《边城》
42、戏曲：曹禺的主要作品《雷雨》，成为中国话剧艺术成熟的标志。
43、国统区: 巴金《第四病室》和《寒夜》
沦陷区:张爱玲《金锁记》
抗日根据地:赵树理《小二黑结婚》 丁玲《太阳照在桑干河上》 孙犁《荷花淀》
44、戏剧 《屈原》是郭沫若这一时期历史剧当中成就最高、影响最大的代表作。
45、杜鹏程《保卫延安》是当代第一部战争小说。
46、伤痕小说 以刘心武的《班主任》为发端，以卢新华的《伤痕》为标志。
47、先锋小说 《冈底斯的诱惑》、《山上的小屋》
48、寻根文学 阿城 《棋王》、《孩子王》、《树王》 王安忆《小鲍庄》
张承志《黑骏马》、《北方的河》 韩少功《爸爸爸》、《女女女》
49、诗歌 “朦胧诗派”以舒婷、顾城、北岛为代表
舒婷《致橡树》 北岛《回答》
50、巴金的《随想录》可以算作20世纪80年代散文的最高成就。
51、新写实小说与新历史小说 以池莉的《烦恼人生》为代表 莫言的《红高粱》
52、荷马史诗包括《伊利昂纪》和《奥德修纪》两部作品。
53、希腊三大悲剧家：埃斯库罗斯“悲剧之父”《被缚的普罗米修斯》
索福克勒斯“戏剧艺术的荷马”《俄狄浦斯王》
欧里庇得斯“希腊悲剧的典范” 《美狄亚》
54、但丁是中世纪欧洲文坛上最伟大的作家，意大利民族文学的奠基者。 代表作《神曲》集
中世纪申雪、哲学、科学和文化艺术之大成，是第一部用意大利语写成的作品。
55、意大利 薄伽丘《十日谈》 拉伯雷《巨人传》 塞万提斯《堂·吉诃德》
56、乔叟 被誉为“英国诗歌之父”。代表作《坎特伯雷故事集》。
57、 莎士比亚 被誉为“英国戏剧之父” 《罗密欧与朱丽叶》
四大悲剧：《哈姆莱特》 《奥赛罗》 《李尔王》 《麦克白》
58、法国 莫里哀 《吝啬鬼》、《伪君子》
59、德国 歌德 《少年维特之烦恼》、《浮士德》
席勒 《阴谋与爱情》
60、古代亚非文学：这一阶段的亚非文学代表作品有古希伯来文学的总汇《旧约》、古代印度诗
歌总集《吠陀》几史诗《摩诃婆罗多》和《罗摩衍那》。
61、中古亚非文学：中国、印度、阿拉伯、日本、朝鲜等国的文学皆取得了极其辉煌的成就。
阿拉伯文学的《一千零一夜》
62、雨果 法国浪漫主义文学运动的领袖。他以“美丑对照”的原则创作《巴黎圣母院》 《悲惨
世界》是后期创作的代表作。
63、普希金：既是俄国浪漫主义文学的主要代表，又是俄国现实文学的奠基人，被誉为“俄国文
学之父”、“伟大的俄国人民诗人”。
长篇诗体小说《叶甫盖尼 奥涅金》是俄国现实主义文学的奠基之作。
64、法国浪漫主义作家：大仲马《三个火枪手》、司汤达《红与黑》、巴尔扎克《人间喜剧》、《高
老头》 福楼拜《包法利夫人》
65、果戈里是俄国“自然派”的创始人，俄国批判现实主义文学的杰出代表。
讽刺喜剧《钦差大臣》 和长篇小说《死魂灵》
66、屠格涅夫是19世纪俄国杰出的现实主义作家，与列夫·托尔斯泰、陀思妥耶夫斯基在世界
文坛被誉为“俄国文学三杰”。
67、易卜生：挪威著名戏剧家，被誉为“现代戏剧之父”。 《玩偶之家》
68、海明威 《老人与海》集中体现了海明威文学创作的最高成就。
69、现代主义文学主要作家有艾略特《荒原》、卡夫卡、奥尼尔、乔伊斯《尤利西斯》和福克纳
《喧哗与骚动》。
《喧哗与骚动》和法国作家普斯特的《追忆逝水年华》、爱尔兰小说家乔伊斯的《尤利西斯》并
称为意识流小说的三大杰作。
70、应用文的含义：是应用写作的表现形态，是写作学的一个重要的分支，是机关、团体、企
事业单位和个人在日常工作和学习中，班里各种十五而形成并使用的、具有使
用价值和一定惯用文章体式的文字信息载体。
71、应用文的特点：实用性 真实性 针对性 政策性和强制性 时效性 具有固定的格式和
处理程序
72、应用文语言的特点：准确 简练 质朴 规范
73、应用文的表达方式：记叙 说明 议论
74、计划的概念 计划是单位或个人为完成工作、学习等任务，预先制定出任务具体内容、措施
和步骤的应用文体。
75、计划的种类（按内容划分）：综合性计划 专项计划
76、通报：用来表彰先进，批评错误，传达重要指示精神或情况时使用的公文文种。
77、通知：指需要特定机关和人员知道、办理事宜时使用的一种公文文种。
78、请示 下级机关请求上级机关领导、业务指导只是一件和批准事项是使用的一种陈述性文件。
79、消息：即消息报道，也叫新闻，是以简洁明快的语言，对现实生活当中发生的具有传播价
值的事实，进行迅速及时报道的文体。
80、消息的要素：“五个W”，即when（时间）、where（地点）、who（人物）、what（事情）、why
（为什么）
81、消息的结构：1.倒金字塔式 2.金字塔式结构 3.倒金字塔与金字塔结合的结构
篇三：北航数值分析大作业3
一、算法设计方案
1.使用牛顿迭代法，对原题中给出的xi?0.08i，yj?0.5?0.05j，
,0?j?20)的11*21组xi,yj分别求出原题中方程组的一组解，于是得(0?i?10
到一组和xi,yi对应的ui,tj。
2.对于已求出的ui,tj，使用分片二次代数插值法对原题中关于z,t,u的数表进行插值得到
zij。于是产生了z=f(x,y)的11*21个数值解。
3.从k=1开始逐渐增大k的值，并使用最小二乘法曲面拟合法对z=f(x,y)进行拟合，得到每次的k,?。当??10时结束计算，输出拟合结果。
**
4.计算f(xi*,y*),p(x,y,2,???,8,j?1,2,???,5)的值并输出结果，以观察p(x,y)逼jij)(i?1**
近f(x,y)的效果。其中xi?0.1i,yj?0.5?0.2j。
?7
二、算法实现方案 1、求f(x,y)：
（1）Newton法解非线性方程组
?0.5cost?u?v?w?x?2.67?t?0.5sinu?v?w?y?1.07?
?????(1)， ?
?0.5t?u?cosv?w?x?3.74??t?0.5u?v?sinw?y?0.79
其中，t, u, v ,w为待求的未知量，x, y为代入的已知量。
设??(t,u,v,w)T，给定精度水平?1?10?12和最大迭代次数M，则解该线性方程组的迭代格式为：
?在?*附近选取初值?(0)?(t(0),u(0),v(0),w(0))T，?(k+1)
??(k)?F?(?(k))?1F(?(k))??
?k?
0,1,?
(k)(k?1)
迭代终止条件为???
?
/?(k)
?
??1，若k?M时仍未达到迭代精度，则迭代计算失
败。
其中，雅可比矩阵
?0.5*cos(t) + u + v + w - x - 2.67??t + 0.5*sin(u) + v + w - y - 1.07?
?， F(?)??
?0.5*t + u + cos(v) + w - x - 3.74????t + 0.5*u + v + sin(w) - y - 0.79?
111???0.5*sin(t)
??10.5*cos(u)11?， F?(?)???0.51?sin(v)1???
10.51cos(w)??
（2）分片双二次插值：
根据题目给出的表格，
ti?ih??(i?0,1,ui?j??(j?0,1,
对于给定的(t,u),如果(t,u)满足
,5),5)
???????） （其中，h?0.2，
ti?
hh
?t?ti?,?????????i?322
uj?
?
?u?uj?,??????j?322
?
则选择(tk,ur)??(k?i?1,i,i?1;r?j?1,j,j?1)为插值节点，相应的插值多项式为
h(t,u)?
其中，
k?i?1r?j?1
??l(t)l(u)g(t,u)??????????
k
r
k
r
i?1j?1
lk(t)?lr(u)?
t?tm
(k?i?1,i,i?1)?m?i?1tk?tm
m?k
i?1
u?un
????(r?j?1,j,j?1)?n?j?1ur?un
n?r
j?1
如果t?t1?
hh??
或t?t4?,则在式（2）中取i=1或i=4; 如果u?u1?或u?u4?,2222
则在式（2）中取u=1或u=4。
在区域D?{(x,y)|0?x?0.8,0.5?y?1.5}上，将
(xi,yj)?(xi?0.08i,i=0,1,…,10;yj?0.5+0.2j,j=0,1,…20)代入到非线性方程组（1）中，用
Newton法解出(ti,j,ui,j)，再由分片双二次插值得zi,j?h(ti,j,ui,j)，则有zi,j?f(xi,yj)（i=0,1,…,10;j=0,1,…,20），即求出了z?f(x,y)。 2、求p(x,y)：
乘积型最小二乘拟合曲面：
（1）求系数矩阵C：
C?(BTB)?1BTUG(GTG)?1
其中，
?1x0?
1x1
B??
??
?1x10??1y0?
1y1
G??
????1y20
x02
x12x102y02y12x202
x0k??x1k?
??x10k??y0k?
?y1k?
??y20k??
?z0,0
?z1,0
U??
???z10,0
z0,1z1,1z10,1
z0,20?z1,20??,z?f(x,y)??(i?0,1,
ij
?i,j?z10,20?
,10;j?0,1,,20)
计算中涉及到对矩阵求逆，接着在后面将会具体说明列主元的高斯消去法求矩阵的逆的方法。
（2）确定最小的k值，拟合曲面：
p(x,y)?
设?
10
20
r,s?0
?c
k
rsxyrs
???(f(xi,yj)?p(xi,yj))2，给定精度水平?2?10?7和最大迭代次数N，
i?0j?0
则确定最小k值的迭代格式为：
k
?(k)rs
?p(xi,yj)??
crsxiyj???(i?0,1,,10;j?0,1,
r,s?0
?
?(k)1020
(k)2
?????(f(xi,yj)?p(xi,yj))
i?0j?0?
?k?0,1,2,??
,20)
迭代终止条件为?
(k)
??2，若k?N时仍未达到迭代精度，则迭代计算失败。
待确定满足精度条件的最小k值后，就可以进行曲面拟合计算了。 3、关于列主元的高斯消去法求矩阵的逆：
设非奇异矩阵A?Cn?n，且B?A ，则
?1
AB?I，
对B和I列分块，有
A(b1,
即
,bj,,bn)?(e1,,ej,,en)
Abj?ej,??j?1,2,
其中，
,n
ej?(0,
j
,0,1,0,
,0)T
应用列主元的高斯消去法线性解方程组 则B?(b1,
Abj?ej,??j?1,2,
,n，
,bj,,bn)即为A的逆。
注：若A不可逆，则此算法失效。 三、源程序
#include "stdafx.h" #include
const int M= 500;//迭代最大次数
const double E=1.0e-12;//确定牛顿迭代精度水平 const double E1=1.0e-7;//确定拟合精度水平 const int kmax=9;//k的最大值 const double matrix[6][6]={ };
const double mat_t[6]={0, 0.2, 0.4, 0.6, 0.8, 1.0}; const double mat_u[6]={0, 0.4, 0.8, 1.2, 1.6, 2.0}; double U[11][21];//对f(x,y)的值进行存储 double C[kmax+1][kmax+1];//拟合系数矩阵
////////////两数绝对值取大///////////// double max(double x,double y){ }
/////////高斯消元法解线性方程组///////
void linear_solution(double f[4],double ff[4][4], double (&delta)[4]) {
return fabs(x)>fabs(y)?fabs(x):fabs(y); {-0.5, -0.34, 0.14, 0.94, 2.06, 3.5}, {-0.42, -0.5, -0.26, 0.3, 1.18, 2.38}, {-0.18, -0.5, -0.5, -0.18, 0.46, 1.42}, {0.22, -0.34, -0.58, -0.5, -0.1, 0.62}, {0.78, -0.02, -0.5, -0.66, -0.5, -0.02}, {1.5, 0.46, -0.26, -0.66, -0.74, -0.5}
}
int i, j, k, ik; double tmp, mik; for(k=0; k<3; k++) { }
delta[3]=f[3]/ff[3][3]; for(k=2;k>=0;k--) { }
tmp=0;
tmp+=ff[k][j]*delta[j]; ik=k;
for(i=k;i<4;i++) {
tmp=ff[k][j]; ff[k][j]=ff[ik][j]; ff[ik][j]=tmp; }
tmp=f[k]; f[k]=f[ik]; f[ik]=tmp;
for(i=k+1;i<4;i++) { }
mik=ff[i][k]/ff[k][k]; for(j=k;j<4;j++)
ff[i][j]=ff[i][j]-mik*ff[k][j]; if(fabs(ff[ik][k]) ik=i; for(j=k;j<4;j++) f[i]-=mik*f[k]; for(j=k+1;j<4;j++) delta[k]=(f[k]-tmp)/ff[k][k]; ////////Newton法求解非线性方程组//////// void equation_solution(double x, double y,double &u,double &t){ double v=1.0,w=1.0; u=1.0; t=1.0; double delta[4],f[4],ff[4][4]; int k=0; for(k=0;k<=M;k++) 篇四：北航数值分析大作业3 《数值分析B》 第三次数值分析大作业 院系：04 能源与动力工程学院 姓名： 王 开 逍 学号： SY1104207 一．算法设计方案： 1、使用牛顿迭代法，对原题中给出的X(I)=0.08*I，Y(J)=0.5+0.05*J的11*21组X(I),Y(J)分别求得原题非线性方程组的一组解，于是得到一对和X(I),Y(J)对应的T(I),U(J). 2、对于已经求出来的U(I)，T(J)，使用分片二次代数插值法对原题中关于Z,T,U的数表进行插值，得到Z(I,J).于是产生了Z=F(X,Y)的11*21个数值解. 3、从K=1开始逐渐增大K的值，并使用最小二乘曲面拟合法对Z=F(X,Y)进行拟合，得到每次的K和精度TAO，当TAO 4、计算F(X*,Y*),P(X*,Y*)的值，以观察P(X,Y)逼近F(X,Y)的效果，其中X*(I)=0.1*I,Y*(J)=0.5+0.2*J. 二．源程序： 该计算程序使用FORTRAN语言编写 module linear3 !创建模块LINEAR3,用于被主程序调用 use imsl implicit none contains real(kind=8) function max12(a,b,c,d) !求向量的无穷范数的子程序 real(kind=8) a,b,c,d real(kind=8) a0,b0,c0,d0 a0=abs(a) b0=abs(b) c0=abs(c) d0=abs(d) if(a0>b0)then b0=a0 end if if(b0>c0)then c0=b0 end if if(c0>d0)then d0=c0 end if max12=d0 end function subroutine qiugeng(x,y,t,u) !求题中所给非线性方程组的解，并且使用牛顿迭代法 real(kind=8) x,y,t,u real(kind=8) w,v,a(4,4),det(4),b(4) integer i,j t=0.5 u=0.5 w=1 v=1 do while(1>0) a(1,1)=0-0.5*sin(t) a(2,2)=0.5*cos(u) a(3,3)=0-sin(v) a(4,4)=cos(w) do i=1,4 do j=1,4 if(i/=j)then a(i,j)=1 end if end do end do a(3,1)=0.5 a(4,2)=0.5 b(1)=2.67-0.5*cos(t)-u-v-w+x b(2)=1.07-t-0.5*sin(u)-v-w+y b(3)=3.74-0.5*t-u-cos(v)-w+x b(4)=0.79-t-0.5*u-v-sin(w)+y det=a.ix.b if(max12(det(1),det(2),det(3),det(4))/max12(t,u,v,w)<=0.1e-11) exit t=t+det(1) u=u+det(2) v=v+det(3) w=w+det(4) end do end subroutine qiugeng real(kind=8) function p22(x,y,x0,y0,h,tao,z) !使用分片二次插值法，对Z=F(X,Y)进行插值 real(kind=8) x,y,x0(:),y0(:),h,tao,z(:,:) real(kind=8) lx(3),ly(3) integer i,j,n,m,k1,k2 m=size(y0,1) n=size(x0,1) do i=3,n-2 if(x>(x0(i)-0.5*h) .and. x<=(x0(i)+0.5*h)) then k1=i end if end do if(x<=(x0(2)+0.5*h)) then k1=2 end if if(x>(x0(n-1)-0.5*h)) then k1=n-1 end if do i=3,m-2 if(y>(y0(i)-0.5*tao) .and. y<=(y0(i)+0.5*tao)) then k2=i end if end do if(y<=(y0(2)+0.5*tao)) then k2=2 end if if(y>(y0(m-1)-0.5*tao)) then k2=m-1 end if lx(1)=(x-x0(k1))*(x-x0(k1+1))/((x0(k1-1)-x0(k1))*(x0(k1-1)-x0(k1+1))) lx(2)=(x-x0(k1-1))*(x-x0(k1+1))/((x0(k1)-x0(k1-1))*(x0(k1)-x0(k1+1))) lx(3)=(x-x0(k1))*(x-x0(k1-1))/((x0(k1+1)-x0(k1))*(x0(k1+1)-x0(k1-1))) ly(1)=(y-y0(k2))*(y-y0(k2+1))/((y0(k2-1)-y0(k2))*(y0(k2-1)-y0(k2+1))) ly(2)=(y-y0(k2-1))*(y-y0(k2+1))/((y0(k2)-y0(k2-1))*(y0(k2)-y0(k2+1))) ly(3)=(y-y0(k2-1))*(y-y0(k2))/((y0(k2+1)-y0(k2))*(y0(k2+1)-y0(k2-1))) p22=0.0 do i=1,3 do j=1,3 p22=p22+lx(i)*ly(j)*z(k2-2+j,k1-2+i) end do end do end function real(kind=8) function pxy(x,y,c) !关于F(X,Y)的X,Y近似多项式表达式的值 real(kind=8) x,y,c(:,:) integer i,j,m,n real(kind=8) t,p pxy=0.0 m=size(c,1) n=size(c,2) do i=1,m do j=1,n t=x**(i-1) p=y**(j-1) pxy=pxy+c(i,j)*t*p end do end do end function end module program main !主程序开始的地方 use linear3 !调用上面的LINEAR3模块 use imsl !调用IMSL库函数，用于曲面拟合的矩阵计算 implicit none real(kind=8),allocatable:: u(:,:),t(:,:),x(:),y(:),u0(:),t0(:),f(:,:) real(kind=8),allocatable:: b(:,:),g(:,:),bt(:,:),gt(:,:),btb(:,:),gtg(:,:),c(:,:) real(kind=8),allocatable:: btuu(:,:),uu(:,:),btbni(:,:),gtgni(:,:),bniu(:,:),ggni(:,:) real(kind=8) z(6,6),xn(8),yn(5) real(kind=8) h,tao,jingdu integer i,j,k allocate(u(11,21),t(11,21),x(11),y(21),u0(6),t0(6),f(11,21)) data z /-0.5,-0.42,-0.18,0.22,0.78,1.5,-0.34,-0.5,-0.5,-0.34,-0.02,0.46,& &0.14,-0.26,-0.5,-0.58,-0.5,-0.26,0.94,0.3,-0.18,-0.5,-0.66,-0.66,& &2.06,1.18,0.46,-0.1,-0.5,-0.74,3.5,2.38,1.42,0.62,-0.02,-0.5/ open(unit=30,file='output/thicknessandf.txt') !将程序运行的结果存储在一个TXT的文档里面 h=0.4 tao=0.2 do i=1,11 x(i)=0.08*(i-1) end do do j=1,21 y(j)=0.5+0.05*(j-1) end do do i=1,11 do j=1,21 call qiugeng(x(i),y(j),t(i,j),u(i,j)) !调用非线性方程组的牛顿迭代解法 end do end do do i=1,6 u0(i)=0+(i-1)*0.4 t0(i)=0+(i-1)*0.2 end do write(30,"('x(i),y(j),f(x,y)的值为：')") do i=1,11 do j=1,21 f(i,j)=p22(u(i,j),t(i,j),u0,t0,h,tao,z) !对于X,Y对应的U,T进行二元插值 write(*,"(f4.2,f6.2, e20.12)") x(i),y(j),f(i,j) write(30,"(f4.2,f6.2, e20.12)") x(i),y(j),f(i,j) end do 篇五：北航数值分析大作业题目三 《数值分析》第三次大作业 一、算法的设计方案： （一）、总体方案设计： （1）解非线性方程组。将给定的(xi,yi)当作已知量代入题目给定的非线性方程组，求得与(xi,yi)相对应的数组t[i][j],u[i][j]。 （2）分片二次代数插值。通过分片二次代数插值运算，得到与数组t[11][21],u[11][21]]对应的数组z[11][21]，得到二元函数z=f(xi,yi)。 （3）曲面拟合。利用x[i],y[j],z[11][21]建立二维函数表，再根据精度的要求选择适当k值，并得到曲面拟合的系数矩阵C[r][s]。 （4）观察和p(xi,yi)的逼近效果。观察逼近效果只需要重复上面（1）和（2）的过程，得到与新的插值节点(xi,yi)对应的f(xi,yi)，再与对应的p(xi,yi)比较即可，这里求解 p(xi,yi)可以直接使用（3）中的C[r][s]和k。 （二）具体算法设计： （1）解非线性方程组 牛顿法解方程组F(x)?0的解x，可采用如下算法： 1）在x附近选取x * (0) * ?D，给定精度水平??0和最大迭代次数M。 2）对于k?0,1,?M执行 ① 计算F(x (k) )和F?(x(k))。 (k) ② 求解关于?x 的线性方程组 (k) F?(x③ 若?x④ 计算x (k) ? )?x(k)??F(x(k)) ??，则取x*?x(k)，并停止计算；否则转④。 x(k) ? (k?1) ?x(k)??x(k)。 ⑤ 若k?M，则继续，否则，输出M次迭代不成功的信息，并停止计算。 （2）分片双二次插值 给定已知数表以及需要插值的节点，进行分片二次插值的算法： 设已知数表中的点为： ? ??xi?x0?ih(i?0,1,?,n) ，需要插值的节点为(x,y)。 y?y?j?(j?0,1,?,m)?0?j 1) 根据(x,y)选择插值节点(xi,yj)： 若x?x1? hh 或x?xn?1?，插值节点对应取i?1或i?n?1， 22 若y?y1? ? 2 或y?yn?1? ? 2 ，插值节点对应取j?1或i?m?1。 hh?x??x?x?,2?i?n?2i??i22 若? ???y??y?y?,2?j?m?2jj??22 则选择(xk,yr)(k?i?1,i,i?1;r?j?1,j,j?1)为插值节点。 2）计算 lk(x)?? t?kj?1 x?xt (k?i?1,i,i?1)x?xt?i?1kt i?1 y?yt l?r(y)??(r?j?1,j,j?1) t?j?1yr?yt t?r 插值多项式的公式为： p(x,y)? k?i?1r?j?1 ??l(x)l?(y)f(x,y) k r k r i?1j?1 注：本步进行插值运算的是(t,u)，利用(xi,yj)与(t,u)的对应关系就可以得到z与 (xi,yj)的对应关系。 （3）曲面拟合 根据插值得到的数表xi,yj,f(xi,yj)进行曲面拟合的过程： 1） 根据拟合节点和基底函数写出矩阵B和G： ?(x0)0?0(x)1 B?? ????(x)0?n?(y0)0(x0)1?(x0)k? ?? (y1)0(x1)1?(x1)k?? G??????? ??1k??(y)0 (xn)?(xn)??m T ?1 T T ?1 (y0)1?(y0)k? ? (y1)1?(y1)k? ?????1k?(ym)?(ym)? 2） 计算 C?(BB)BUG(GG)。 在这里，为了简化计算和编程、避免矩阵求逆，记： A?(BTB)?1BTU，DT?G(GTG)?1 对上面两式进行变形，得到如下两个线性方程组： (BTB)A?BTU，(GTG)D?GT 通过解上述两个线性方程组，则有：C?AD 3） 对于每一个(xi,yj)， T p(xi,yj)???Crs(xi)r(yj)s。 * r?0s?0 kk 4） 拟合需要达到的精度条件为： ?? ??[p(x,y)?u * i j i?0j?0 nm ij ]2?10?7。 其中uij对应着插值得到的数表xi,yj,f(xi,yj)中f(xi,yj)的值。 5） 让k逐步增加，每一次重复执行以上几步，直到 ????[p*(xi,yj)?uij]2?10?7 成立。此时的k值就是要求解最小的k。 i?0j?0 nm 二、源程序： #include #define Epsilon1 1e-12 /*解线性方程组时近似解向量的精度*/ #define M 200 /*解线性方程组时的最大迭代次数*/ #define N 10 /*求解迭代次数时假设的k的最大值，用于定义包含k的存储空间*/ void Newton(); /*牛顿法求解非线性方程组子程序*/ void fpeccz(); /*分片二次代数插值子程序*/ void qmnh(); /*曲面拟合子程序*/ void duibi(); /*对比??和p逼近效果的子程序*/ double x[11],y[21],t[11][21],u[11][21];/*定义全局变量*/ double z[11][21],C[10][10]; double kz; void Newton(double x[11],double y[21])/*牛顿法求解非线性方程组子程序*/ { double X[4],dx[4],F[4],dF[4][4],temp,m,fx,fX; int i,j,k,l,p,ik,n; for(i=0;i<=10;i++) { for(j=0;j<=20;j++) { X[0]=1; /*选取迭代初始向量，四个分别代表t，u，v，w*/ X[1]=1; X[2]=1; X[3]=1; n=0; loop1:{ F[0]=0.5*cos(X[0])+X[1]+X[2]+X[3]-x[i]-2.67; F[1]=X[0]+0.5*sin(X[1])+X[2]+X[3]-y[j]-1.07; F[2]=0.5*X[0]+X[1]+cos(X[2])+X[3]-x[i]-3.74; F[3]=X[0]+0.5*X[1]+X[2]+sin(X[3])-y[j]-0.79; /*求解F（x）*/ dF[0][0]=-0.5*sin(X[0]); /*求解F'（x）*/ dF[0][1]=1; dF[0][2]=1; dF[0][3]=1; dF[1][0]=1; dF[1][1]=0.5*cos(X[1]); dF[1][2]=1; dF[1][3]=1; dF[2][0]=0.5; dF[2][1]=1; dF[2][2]=-sin(X[2]); dF[2][3]=1; dF[3][0]=1; dF[3][1]=0.5; dF[3][2]=1; dF[3][3]=cos(X[3]); /*高斯选主元消去法求解Δx*/ for(k=0;k<3;k++) { ik=k; for(l=k;l<=3;l++) {if(dF[ik][k] temp=0; temp=F[ik]; F[ik]=F[k]; F[k]=temp; for(l=k;l<=3;l++) { temp=0; temp=dF[ik][l]; dF[ik][l]=dF[k][l]; dF[k][l]=temp; } for(l=k+1;l<=3;l++) { m=dF[l][k]/dF[k][k]; F[l]=F[l]-m*F[k]; for(p=k+1;p<=3;p++) {dF[l][p]=dF[l][p]-m*dF[k][p];} } /*消去过程*/ } dx[3]=-F[3]/dF[3][3]; for(k=2;k>=0;k--) { temp=0; for(l=k+1;l<=3;l++) {temp=temp+dF[k][l]*dx[l]/dF[k][k];} dx[k]=-F[k]/dF[k][k]-temp; } temp=0; for(l=0;l<=3;l++) /*求解矩阵范数，用无穷范数*/ { if(temp fx=temp; temp=0; for(l=0;l<=3;l++) { if(temp fX=temp;