傳送門

由於沒有考慮n<=1的情況T了很久啊。 這題很有意思啊。 考試的時候根本不會，騙了30分走人。 實際上變一個形就可以了。 推導過程有點繁雜。 直接粘題解上的請諒解。 不得不說這個推導很妙。 然後就可以矩陣快速冪優化了。 程式碼：

#include<bits/stdc++.h>
#define ll long long
using namespace std;
ll a,b,n,mod;
struct Matrix{
	ll a[3][3];
	Matrix(){a[0][0]=a[0][1]=a[0][2]=a[1][0]=a[1][1]=a[1][2]=a[2][0]
 
=a[2][1]=a[2][2]=0;}
	inline void init(){a[0][0]=0,a[0][1]=1,a[0][2]=0,a[1][0]=1,a[1][1]=1,a[1][2]=1,a[2][0]=0,a[2][1]=2,a[2][2]=1;}
	inline Matrix operator*(Matrix b){
		Matrix c;
		c.a[0][0]=(a[0][0]*b.a[0][0]%mod+a[0][1]*b.a[1][0]%mod+a[0][2]*b.a[2][0]%mod)%mod;
		c.a[0][1]=(a[0][0]*b.a[0][1]%mod+a[0][1]*b.a[1][1]%mod+a[
 
0][2]*b.a[2][1]%mod)%mod;
		c.a[0][2]=(a[0][0]*b.a[0][2]%mod+a[0][1]*b.a[1][2]%mod+a[0][2]*b.a[2][2]%mod)%mod;
		c.a[1][0]=(a[1][0]*b.a[0][0]%mod+a[1][1]*b.a[1][0]%mod+a[1][2]*b.a[2][0]%mod)%mod;
		c.a[1][1]=(a[1][0]*b.a[0][1]%mod+a[1][1]*b.a[1][1]%mod+a[1][2]*b.a[2][1]%mod)%mod;
		c.a[1][2]=(a[1][0]*b.a[0][2]%mod+a[
 
1][1]*b.a[1][2]%mod+a[1][2]*b.a[2][2]%mod)%mod;
		c.a[2][0]=(a[2][0]*b.a[0][0]%mod+a[2][1]*b.a[1][0]%mod+a[2][2]*b.a[2][0]%mod)%mod;
		c.a[2][1]=(a[2][0]*b.a[0][1]%mod+a[2][1]*b.a[1][1]%mod+a[2][2]*b.a[2][1]%mod)%mod;
		c.a[2][2]=(a[2][0]*b.a[0][2]%mod+a[2][1]*b.a[1][2]%mod+a[2][2]*b.a[2][2]%mod)%mod;
		return c;
	}
};
inline ll ksm(ll a,ll b,ll c,int p){
	Matrix ret,x;
	ret.init(),x.init();
	while(p>=1){
		if(p&1)ret=ret*x;
		x=x*x,p>>=1;
	}
	return (a*ret.a[0][1]%mod+b*ret.a[1][1]%mod+c*ret.a[2][1]%mod)%mod;
}
int main(){
	cin>>a>>b>>mod>>n;
	if(n<=1)puts("1");
	else cout<<ksm(a,b,sqrt((ll)(3+a*b)),n-2);
	return 0;
}